Télécharger la plaquette Cast3M


fr en

Télécharger flacr3.eso

Retour à la liste

Numérotation des lignes : 

flacr3

  1. C FLACR3    SOURCE    CB215821  20/11/25    13:29:00     10792          
  2.       SUBROUTINE FLACR3(EPSCSI,EPSI,DELTAT,MCEN,MFACEL,IRC,IYC,IYINIT,
  3.      &     IYFINA,IVCAR,IDX,MLRMAS,MLRH0K,MLRECO,ICHRET,ICHRY)
  4. C************************************************************************
  5. C
  6. C PROJET            :  CASTEM 2000
  7. C
  8. C NOM               :  FLACR3
  9. C
  10. C DESCRIPTION       :  CREBCOM: modele non-homogene
  11. C                      voir FLACR2
  12. C
  13. C LANGAGE           :  FORTRAN 77 + ESOPE 2000 (avec extensions CISI)
  14. C
  15. C AUTEUR            :  A. BECCANTINI, DM2S/SFME/LTMF
  16. C
  17. C************************************************************************
  18. C
  19. C INPUTS
  20. C
  21. C EPSCSI   : (REAL*8), parametre pour controler si CSI=DY/(YF - YI)
  22. C            \in(0,1). Si  CSI<EPSCSI ou CSI>1+EPSCSI un message de
  23. C            warnin est donné.
  24. C
  25. C EPSI     : parametre du critère CREBCOM (REAL*8)
  26. C
  27. C DELTAT   : pas de temps (REAL*8)
  28. C
  29. C MCEN    : pointeur du MELEME contenant les centres des ELTs
  30. C
  31. C MFACEL  : pointeur du MELEME contenant la correspondence
  32. C           CENTRE-FACE-CENTRE
  33. C
  34. C IRC     : pointeur du CHPOINT  contenant la masse volumique.
  35. C
  36. C IYC     : pointeur du CHPOINT contenant les fractions massiques
  37. C
  38. C IYINIT  : pointeur du CHPOINT contenant la fraction massique initiale
  39. C           de la premiere composante de IYC;
  40. C
  41. C IYFINA  : pointeur du CHPOINT contenant la fraction massique finale de
  42. C           la premiere composante de IYC ;
  43. C
  44. C IVCAR   : pointeur du CHPOINT contenant la vitesse caractéristique de
  45. C           la flamme
  46. C
  47. C IDX     : pointeur du CHPOINT contenant la diemnsion de la maille
  48. C
  49. C MLRMAS  : pointeur du LISTREEL contenant les masses molaires
  50. C
  51. C MLRH0K  : pointeur du LISTREEL contenant les energies des formation à
  52. C           0K
  53. C
  54. C MLRCOE  : pointeur du LISTREEL contenant les coeff. stoch. de la
  55. C           reaction chimique
  56. C
  57. C OUTPUTS :
  58. C
  59. C ICHRET  : pointeur de l'increment de l'energie totale par unité de
  60. C           volume
  61. C
  62. C ICHRY   : pointeur de l'increment des densités massiques
  63. C
  64. C************************************************************************
  65. C
  66. C HISTORIQUE (Anomalies et modifications éventuelles)
  67. C
  68. C HISTORIQUE :
  69. C
  70. C
  71. C************************************************************************
  72. C
  73.       IMPLICIT INTEGER(I-N)
  74.  
  75. -INC PPARAM
  76. -INC CCOPTIO
  77. -INC SMLREEL
  78.       POINTEUR MLRECO.MLREEL, MLRMAS.MLREEL, MLRH0K.MLREEL
  79. -INC SMELEME
  80.       POINTEUR MELCEN.MELEME, MELEFL.MELEME
  81. -INC SMCHPOI
  82.       INTEGER N, NC
  83.       POINTEUR MPOVRO.MPOVAL, MPOVAY.MPOVAL,MPOVYI.MPOVAL,
  84.      &     MPOVYF.MPOVAL, MPOVC.MPOVAL, MPODRE.MPOVAL, MPODRY.MPOVAL,
  85.      &     MPOCSI.MPOVAL,MPOCRI.MPOVAL,MPOVDX.MPOVAL
  86. -INC SMLENTI
  87.       POINTEUR MLECEN.MLENTI
  88. C
  89. C
  90. C**** Variables de COOPTIO
  91. C
  92. C      INTEGER IPLLB, IERPER, IERMAX, IERR, INTERR
  93. C     &        ,IOTER,   IOLEC,  IOIMP,     IOCAR,  IOACQ
  94. C     &        ,IOPER,   IOSGB,  IOGRA,     IOSAU,  IORES
  95. C     &        ,IECHO,   IIMPI,  IOSPI
  96. C     &        ,IDIM
  97. C     &        ,MCOORD
  98. C     &        ,IFOMOD, NIFOUR, IFOUR, NSDPGE, IONIVE
  99. C     &        ,NGMAXY, IZROSF, ISOTYP, IOSCR,LTEXLU
  100. C     &        ,NORINC,NORVAL,NORIND,NORVAD
  101. C     &        ,NUCROU, IPSAUV, IFICLE, IPREFI
  102. C
  103. C**** Les variables
  104. C
  105.       INTEGER MCEN,MFACEL,IRC,IYC,IYINIT,IYFINA,IVCAR,ICHRET,ICHRY
  106.      $     ,IGEOM,ICEN,NCEN,NFAC,NLCF,NGCEG,NGCED,NLCEG,NLCED
  107.      &     ,IESP,NESP,IDX
  108.       REAL*8 EPSI, DELTAT, Y1F, Y1I, Y1, VCSIG, VCSID, VCSI2G, VCSI2D
  109.      &     ,EPS12, DY, DYMAX, DCSI, DCSI1, DY1, RHO, DY2, YMAX, DYMAX1
  110.      &     ,ERRTOL, EPSCSI, CSIG
  111.       CHARACTER*8 TYPE
  112.       PARAMETER(ERRTOL=1.0D-6)
  113. C
  114.       NESP=MLRH0K.PROG(/1)
  115.       MELCEN=MCEN
  116.       MELEFL=MFACEL
  117. C
  118. C**** KRIPAD pour la correspondance global/local de centre
  119. C
  120.       CALL KRIPAD(MELCEN,MLECEN)
  121.       IF(IERR .NE. 0)GOTO 9999
  122. C     SEGINI MLECEN
  123.       SEGACT MELCEN
  124.       NCEN=MELCEN.NUM(/2)
  125.       SEGACT MELEFL
  126.       NFAC=MELEFL.NUM(/2)
  127. C
  128. C**** Lectures de CHPOINT
  129. C
  130.       CALL LICHT(IRC,MPOVRO,TYPE,IGEOM)
  131. C     SEGACT MPOVRO*MOD
  132.       CALL LICHT(IYC,MPOVAY,TYPE,IGEOM)
  133. C     SEGACT MPOVAY*MOD
  134.       CALL LICHT(IYINIT,MPOVYI,TYPE,IGEOM)
  135. C     SEGACT MPOVYI*MOD
  136.       CALL LICHT(IYFINA,MPOVYF,TYPE,IGEOM)
  137. C     SEGACT MPOVYF*MOD
  138.       CALL LICHT(IVCAR,MPOVC,TYPE,IGEOM)
  139. C     SEGACT MPOVYF*MOD
  140.       CALL LICHT(IDX,MPOVDX,TYPE,IGEOM)
  141. C     SEGACT MPOVDX*MOD
  142.       CALL LICHT(ICHRET,MPODRE,TYPE,IGEOM)
  143. C     SEGACT MPODRE*MOD
  144.       CALL LICHT(ICHRY,MPODRY,TYPE,IGEOM)
  145. C     SEGACT MPODRY*MOD
  146. C
  147. C**** Creation du MPOVAL qui contient csi
  148. C     Creation du MPOVAL du critere
  149. C
  150.       N=NCEN
  151.       NC=1
  152.       SEGINI MPOCSI
  153.       SEGINI MPOCRI
  154. C
  155. C**** Calcul de DYMAX pour la premiere espece
  156. C     Controle de la densité positive
  157. C     Controle des fractions massiques (Yi>0, sum_i Y_i < 1)
  158. C
  159.       DYMAX = 0.0D0
  160.       YMAX=0.0D0
  161.       DO ICEN=1,NCEN,1
  162.          Y1F = MPOVYF.VPOCHA(ICEN,1)
  163.          Y1I = MPOVYI.VPOCHA(ICEN,1)
  164.          DY = Y1F - Y1I
  165.          DYMAX=MAX(ABS(DY),DYMAX)
  166.          Y1=MPOVAY.VPOCHA(ICEN,1)
  167.          IF(Y1 .LT. 0.0D0)THEN
  168.             WRITE(IOIMP,*) 'CHPO2 < 0 ???'
  169.             CALL ERREUR(21)
  170.             GOTO 9999
  171.          ELSE
  172.             YMAX=Y1
  173.          ENDIF
  174.          RHO=MPOVRO.VPOCHA(ICEN,1)
  175.          IF(RHO .LT. 0)THEN
  176.             WRITE(IOIMP,*) 'CHPO1 < 0 ???'
  177.             CALL ERREUR(21)
  178.             GOTO 9999
  179.          ENDIF
  180.          RHO=MPOVYI.VPOCHA(ICEN,1)
  181.          IF((RHO .LT. 0) .OR. (RHO .GT. 1))THEN
  182.             WRITE(IOIMP,*) 'CHPO3 < 0 ou CHPO3 > 1???'
  183.             CALL ERREUR(21)
  184.             GOTO 9999
  185.          ENDIF
  186.          RHO=MPOVYF.VPOCHA(ICEN,1)
  187.          IF((RHO .LT. 0.0D0) .OR. (RHO .GT. 1.0D0))THEN
  188.             WRITE(IOIMP,*) 'CHPO4 < 0 ou CHPO4 > 1???'
  189.             CALL ERREUR(21)
  190.             GOTO 9999
  191.          ENDIF
  192.          RHO=MPOVC.VPOCHA(ICEN,1)
  193.          IF((RHO .LT. 0))THEN
  194.             WRITE(IOIMP,*) 'CHPO5 < 0'
  195.             CALL ERREUR(21)
  196.             GOTO 9999
  197.          ENDIF
  198.          DO IESP=2,NESP,1
  199.             Y1=MPOVAY.VPOCHA(ICEN,IESP)
  200.             IF(Y1 .LT. 0.0D0)THEN
  201.                WRITE(IOIMP,*) 'CHPO2 < 0 ???'
  202.                CALL ERREUR(21)
  203.                GOTO 9999
  204.             ELSE
  205.                YMAX=YMAX+Y1
  206.             ENDIF
  207.          ENDDO
  208.          IF(YMAX .GT. 1.0D0)THEN
  209.             WRITE(IOIMP,*) 'sum CHPO2 > 1 ???'
  210.             CALL ERREUR(21)
  211.             GOTO 9999
  212.          ENDIF
  213.       ENDDO
  214. C
  215. C**** Calcul de CSI \in (0,1)
  216. C
  217.       DO ICEN=1,NCEN,1
  218.          Y1F = MPOVYF.VPOCHA(ICEN,1)
  219.          Y1I = MPOVYI.VPOCHA(ICEN,1)
  220.          DY = Y1F - Y1I
  221.          IF(ABS(DY) .LE. (ERRTOL*DYMAX))THEN
  222. C
  223. C********** Pas de combustion en cette region
  224. C
  225.             MPOCSI.VPOCHA(ICEN,1)=0.0D0
  226. C
  227.          ELSE
  228.             Y1 = MPOVAY.VPOCHA(ICEN,1)
  229.             VCSIG = (Y1 - Y1I) / DY
  230. C
  231. C********** On n'accepte pas csi > 1.1 or csi < -0.1
  232. C
  233.             IF((VCSIG .GT. (1.0D0+EPSCSI)) .OR.
  234.      &           (VCSIG .LT. (-1*EPSCSI)))THEN
  235.                WRITE(IOIMP,*) 'Progress variable = ???'
  236. C              21 2
  237. C              Données incompatibles
  238.                CALL ERREUR(21)
  239.                GOTO 9999
  240.             ELSEIF(VCSIG .GT. 1.0D0)THEN
  241.                VCSIG = 1.0D0
  242.             ELSEIF(VCSIG .LT. 0.0D0)THEN
  243.                VCSIG = 0.0D0
  244.             ENDIF
  245.             MPOCSI.VPOCHA(ICEN,1)= VCSIG
  246.          ENDIF
  247.       ENDDO
  248. C
  249. C**** Le critere CREBCOM
  250. C
  251.       DO NLCF = 1, NFAC, 1
  252. C
  253. C******* NLCF  = numero local du centre de facel
  254. C        NGCEG = numero global du centre ELT "gauche"
  255. C        NLCEG = numero local du centre ELT "gauche"
  256. C        NGCED = numero global du centre ELT "droite"
  257. C        NLCED = numero local du centre ELT "droite"
  258. C
  259.          NGCEG = MELEFL.NUM(1,NLCF)
  260.          NGCED = MELEFL.NUM(3,NLCF)
  261.          NLCEG = MLECEN.LECT(NGCEG)
  262.          NLCED = MLECEN.LECT(NGCED)
  263. C
  264.          VCSIG=MPOCSI.VPOCHA(NLCEG,1)
  265.          VCSID=MPOCSI.VPOCHA(NLCED,1)
  266.          VCSI2G=VCSIG*VCSIG
  267.          VCSI2D=VCSID*VCSID
  268. C
  269.          IF(NLCEG .EQ. NLCED)THEN
  270. C
  271. C********** Murs
  272. C
  273.             MPOCRI.VPOCHA(NLCEG,1)=MPOCRI.VPOCHA(NLCEG,1) + (0.5D0 *
  274.      &           VCSI2D)
  275. C
  276.          ELSE
  277. C
  278.             MPOCRI.VPOCHA(NLCEG,1)=MPOCRI.VPOCHA(NLCEG,1) +
  279.      &           (VCSI2D - (0.5D0 * VCSI2G))
  280.             MPOCRI.VPOCHA(NLCED,1)=MPOCRI.VPOCHA(NLCED,1) +
  281.      &           (VCSI2G - (0.5D0 * VCSI2D))
  282. C
  283.          ENDIF
  284.       ENDDO
  285. C
  286. C**** Calcul des increments de l'energie et des densités
  287. C
  288.       EPS12 = EPSI * EPSI
  289.       DO ICEN = 1, NCEN, 1
  290.          CSIG = MPOCSI.VPOCHA(ICEN,1)
  291.          VCSIG =  MPOCRI.VPOCHA(ICEN,1)
  292. C
  293. C******* In 2D, contribution of the ideal upper and lower cells
  294. C
  295.          IF(IDIM .EQ. 2) VCSIG = VCSIG + (CSIG * CSIG)
  296.          IF(VCSIG .GE. (EPS12*(1.0D0+ERRTOL)))THEN
  297. C
  298. C********** Il y a combustion
  299. C
  300.             DCSI=MPOVC.VPOCHA(ICEN,1)*DELTAT/MPOVDX.VPOCHA(ICEN,1)
  301.             DCSI1=1.0D0-CSIG
  302.             DCSI=MIN(DCSI,DCSI1)
  303.             DY1=(MPOVYF.VPOCHA(ICEN,1)-MPOVYI.VPOCHA(ICEN,1))*DCSI
  304. C
  305. C********** On force la positivité de MPOVAY.VPOCHA(ICEN,1)
  306. C
  307.             IF(DY1 .GT. 0.0D0)THEN
  308.                DY2 = 1.0D0 - MPOVAY.VPOCHA(ICEN,1)
  309.             ELSE
  310.                DY2 = MPOVAY.VPOCHA(ICEN,1)
  311.             ENDIF
  312.             DY1=SIGN(MIN(ABS(DY2),ABS(DY1)),DY1)
  313. C
  314.             RHO=MPOVRO.VPOCHA(ICEN,1)
  315.             MPODRY.VPOCHA(ICEN,1)=DY1*RHO
  316.             MPODRE.VPOCHA(ICEN,1)=-1.0D0*DY1*RHO*MLRH0K.PROG(1)
  317.             DY1=DY1/(MLRMAS.PROG(1)*MLRECO.PROG(1))
  318.             DO IESP=2,NESP,1
  319.                DY=DY1*(MLRMAS.PROG(IESP)*MLRECO.PROG(IESP))
  320.                DYMAX1=ABS(DYMAX/(MLRMAS.PROG(1)*MLRECO.PROG(1))*
  321.      &              (MLRMAS.PROG(IESP)*MLRECO.PROG(IESP)))
  322.                IF(DY .GT. 0)THEN
  323.                   DY2=1.0D0 - MPOVAY.VPOCHA(ICEN,IESP)
  324.                ELSE
  325.                   DY2=MPOVAY.VPOCHA(ICEN,IESP)
  326.                ENDIF
  327. C
  328. C************* On force la positivité de MPOVAY.VPOCHA(ICEN,IESP)
  329. C
  330.                IF((ABS(DY)-ABS(DY2)).GE.(0.5D0*DYMAX1))THEN
  331.                   WRITE(IOIMP,*) 'CHPO2, CHPO3, CHPO4 = ???'
  332. C                 Données incompatibles
  333.                   CALL ERREUR(21)
  334.                   GOTO 9999
  335.                ELSE
  336.                   DY=SIGN(MIN(ABS(DY2),ABS(DY)),DY)
  337.                ENDIF
  338.                MPODRY.VPOCHA(ICEN,IESP)=DY*RHO
  339.                MPODRE.VPOCHA(ICEN,1)=MPODRE.VPOCHA(ICEN,1)-(DY*RHO
  340.      $              *MLRH0K.PROG(IESP))
  341.             ENDDO
  342.          ENDIF
  343.       ENDDO
  344. C
  345.       SEGDES MELCEN
  346.       SEGDES MELEFL
  347.       SEGSUP MLECEN
  348. C
  349.       SEGDES MPOVRO
  350.       SEGDES MPOVAY
  351.       SEGDES MPOVYI
  352.       SEGDES MPOVYF
  353.       SEGDES MPOVC
  354.       SEGDES MPOVDX
  355.       SEGDES MPODRE
  356.       SEGDES MPODRY
  357.       SEGSUP MPOCSI
  358.       SEGSUP MPOCRI
  359. C
  360.  9999 RETURN
  361.       END
  362.  
  363.  
  364.  
  365.  
  366.  
  367.  
  368.  
  369.  
  370.  
  371.  
  372.  
  373.  
  374.  
  375.  
  376.  

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