1 : C FLAM2 SOURCE CHAT 05/01/13 00:03:11 5004 2 : SUBROUTINE FLAM2(N,RHO,CV,R,YH2,YO2,YN2,YH2O,T,ZK,ZE,YH2u, 3 : & ZH2,ZO2,ZN2,ZH2O,Q,CEBU,Dt,CMIN,NCVOPT) 4 : C 5 : C--------------------------------------------------------------------- 6 : C Combustion turbulente de l'hydrogène : Modèle Eddy Break Up (EBU) 7 : C--------------------------------------------------------------------- 8 : C E/ N : Nombre d'états à considérer (dimension des vecteurs) 9 : C E/ RHO : Densité du mélange 10 : C /S CV : Chaleur spécifique à volume constant 11 : C /S R : Constante des gaz du mélange 12 : C E/ YH2, ... Fraction massique de H2, ... 13 : C E/ T : Température du mélange 14 : C E/ ZK : Energie cinétique turbulente 15 : C E/ ZE : Dissipation de l'énergie cinétique turbulente 16 : C E/ YH2u : Fraction massique initiale d'hydrogène 17 : C /S ZH2, ... Fraction massique de H2, ... après combustion 18 : C /S Q : Energie libérée par la combustion 19 : C E/ CEBU : "Constante" du modèle EBU 20 : C E/ Dt : Pas de temps considéré 21 : C E/ CMIN : Seuil relatif pour le controle de la combustion 22 : C E/ NCVOPT : Flag pour calcul cv : 1 si linéaire et 2 si quadratique 23 : C--------------------------------------------------------------------- 24 : C 3/11/99 : Effets de seuil sur la variable de progrès via CMIN 25 : C--------------------------------------------------------------------- 26 : C 27 : IMPLICIT INTEGER(I-N) 28 : IMPLICIT REAL*8 (A-H,O-Z) 29 : C 30 : DIMENSION RHO(*),YH2(*),YO2(*),YN2(*),YH2O(*),YH2u(*) 31 : DIMENSION CV(*),T(*),Q(*),R(*) 32 : DIMENSION ZK(*),ZE(*) 33 : DIMENSION ZH2(*),ZO2(*),ZN2(*),ZH2O(*) 34 : C 35 : REAL*8 MH2,MO2,MN2,MH2O 36 : DATA Rg,MH2,MO2,MN2,MH2O/8.313D0,2.0D-3,32.0D-3,28.0D-3,18.0D-3/ 37 : DATA Cf,Ta/3.3D5,1.0D4/ 38 : DATA RRH2,RRO2,RRN2,RRH2O/4156.5D0,259.8D0,296.9D0,461.8D0/ 39 : DATA ho_H2O/-1.3272D7/ 40 : C 41 : C- Initialisation des fractions massiques et bornage de la turbulence 42 : C 43 : DO I=1,N 44 : ZK(I) = ABS(ZK(I)) + 1.D-20 45 : ZE(I) = ABS(ZE(I)) + 1.D-20 46 : YH2(I) = MAX(YH2(I),0.D0) 47 : YH2(I) = MIN(YH2(I),YH2u(I)) 48 : YO2(I) = ABS(YO2(I)) 49 : YH2O(I) = ABS(YH2O(I)) 50 : YN2(I) = 1.D0 - YH2(I) - YO2(I) - YH2O(I) 51 : ENDDO 52 : C 53 : C--------------------------- 54 : C- Résolution de Dc/Dt = -wc 55 : C--------------------------- 56 : C 57 : DO I=1,N 57 : CALL XCV(T(I),CVH2,CVN2,CVO2,CVH2O,NCVOPT) 58 : CV(I) = YH2(I)*CVH2 + YO2(I)*CVO2 + YN2(I)*CVN2 + YH2O(I)*CVH2O 59 : R(I) = YH2(I)*RRH2 + YO2(I)*RRO2 + YN2(I)*RRN2 + YH2O(I)*RRH2O 60 : * 61 : RH2 = YH2(I)*RHO(I) 62 : RO2 = YO2(I)*RHO(I) 63 : RH2O = YH2O(I)*RHO(I) 64 : ROPREC = RH2 65 : VPROG = 1.D0 - YH2(I)/YH2U(I) 66 : C 67 : C- On empeche le démarrage de la combustion pour des variations de 68 : C- densité de H2 inférieur à CMIN fois la densité de H2 initiale 69 : C 70 : IF (VPROG.LE.CMIN) THEN 71 : Q(I) = 0.D0 72 : ELSE 73 : C 74 : C- Cinétique chimique du modèle : Droh2/Dt = max (Arrhenius,EBU) 75 : C 76 : A_EBU = CEBU*ZE(I)/ZK(I)*VPROG 77 : A_ARR = 2.0D0*Cf*RO2/MO2*EXP(-Ta/T(I)) 78 : A = MAX(A_EBU,A_ARR) 79 : C 80 : C- Si la variable de progrès VPROG vaut 1 à CMIN près, tout brule. 81 : C 82 : IF (VPROG.GE.(1.D0-CMIN)) THEN 83 : RH2 = 0.D0 84 : ELSE 85 : RH2 = ROPREC/(1.0D0 + A*DT) 86 : ENDIF 87 : C 88 : C- Energie libérée et spéciation pour une masse d'H2 brulée égale à DRHO 89 : C- Correction éventuelle de DRHO pour le cas où déficit en O2 90 : C 91 : DRHO = RH2-ROPREC 92 : RO2PRE = RO2 93 : RO2 = RO2 + MO2/(2.0D0*MH2)*DRHO 94 : IF (RO2.LT.0.D0) THEN 95 : RO2 = 0.D0 96 : DRHO = -2.0D0*MH2/MO2*RO2PRE 97 : RH2 = ROPREC + DRHO 98 : ENDIF 99 : RH2O = RH2O - MH2O/MH2*DRHO 100 : Q(I) = MH2O/MH2*DRHO*ho_H2O/DT 101 : ENDIF 102 : C 103 : C- Répartition entre les différentes espèces après combustion 104 : C 105 : ZH2(I) = RH2/RHO(I) 106 : ZO2(I) = RO2/RHO(I) 107 : ZH2O(I) = RH2O/RHO(I) 108 : ZN2(I) = 1.D0 - ZH2(I) - ZO2(I) - ZH2O(I) 109 : ENDDO 110 : C 111 : RETURN 112 : END 113 : 114 : 115 :
© Cast3M 2003 - All rights reserved.
Disclaimer