Télécharger flamarrh.dgibi

Retour à la liste

Numérotation des lignes :

  1. * fichier : flamarrh.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4. **********************************************************************
  5. **** ****
  6. **** OPERATEUR FLAM ****
  7. **** Combustion Hydrogene-Oxygene, cinetique de type Arrhenius ****
  8. **** ****
  9. **** A. BECCANTINI DRN/DMT/SEMT/LTMF FEVRIER 1999 ****
  10. **********************************************************************
  11.  
  12. 'OPTION' 'DIME' 2 ;
  13. 'OPTION' 'ELEM' QUA4 ;
  14. 'OPTION' 'TRAC' 'X' ;
  15. 'OPTION' 'ECHO' 0 ;
  16.  
  17. *
  18. *** GRAPH
  19. *
  20.  
  21. * GRAPH = VRAI ;
  22. GRAPH = FAUX ;
  23.  
  24. ***********************************************************************
  25. ***************** CALCUL DU TAUX DE REACTION *************************
  26. ***********************************************************************
  27.  
  28. 'DEBPROC' TAUX ;
  29. *
  30. **** La meme chose que FLAM ARRHENIUS
  31. *
  32.  
  33. ARGUMENT TS*CHPOINT A*FLOTTANT b*FLOTTANT con*FLOTTANT Ta*FLOTTANT
  34. H0H2*FLOTTANT H0O2*FLOTTANT H0H2O*FLOTTANT DT*FLOTTANT
  35. RYH2*CHPOINT RYO2*CHPOINT RYH2O*CHPOINT T*CHPOINT ;
  36.  
  37.  
  38. MH2 = 2.016D-3 ;
  39. MO2 = 31.998D-3 ;
  40. MH2O = 18.015D-3 ;
  41.  
  42. MAIL = 'EXTRAIRE' TS 'MAILLAGE' ;
  43. NPOI = 'NBEL' MAIL ;
  44.  
  45. LIRYH2 = 'PROG' NPOI * 0.0D0 ;
  46. LIRYO2 = 'PROG' NPOI * 0.0D0 ;
  47. LIRYH2O = 'PROG' NPOI * 0.0D0 ;
  48. LIDE = 'PROG' NPOI * 0.0D0 ;
  49.  
  50. 'REPETER' BL1 NPOI ;
  51. P1 = 'POIN' &BL1 MAIL ;
  52. RYH2I = 'EXTRAIRE' RYH2 'H2 ' P1 ;
  53. RYO2I = 'EXTRAIRE' RYO2 'O2 ' P1 ;
  54. RYH2OI = 'EXTRAIRE'RYH2O 'H2O ' P1 ;
  55. TS0 = 'EXTRAIRE' TS 'SCAL' P1 ;
  56. T0 = 'EXTRAIRE' T 'SCAL' P1 ;
  57. XH2 = RYH2I '/' MH2 ;
  58. XO2 = RYO2I '/' MO2 ;
  59. 'SI' (T0 > TS0) ;
  60. 'SI' (XH2 > (2.0D0 * XO2)) ;
  61. * H2 majoritaire
  62. LAMBDA = EXP(-1.0D0*TA/T0);
  63. LAMBDA = LAMBDA * (T0 ** (-1.0D0*B));
  64. LAMBDA = LAMBDA * (RYH2I ** CON);
  65. LAMBDA = LAMBDA * MO2 * A;
  66. NRYO2I = RYO2I / (1.0D0 + (LAMBDA * DT));
  67. DRYO2 = NRYO2I '-' RYO2I;
  68. DMOL = DRYO2 / MO2;
  69. DRYH2 = 2.0D0 * DMOL * MH2;
  70. DRYH2O = -2.0D0 * DMOL * MH2O;
  71. NRYH2I = RYH2I '+' DRYH2 ;
  72. NRYH2OI = RYH2OI '+' DRYH2O ;
  73. 'MESSAGE' ('CHAINE' 'POINT ' &BL1) ;
  74. 'MESSAGE' 'O2 minoritaire' ;
  75. 'MESSAGE' ('CHAINE' 'DTc = ' (1.0 '/' LAMBDA) ' DT = ' DT) ;
  76. 'SINON' ;
  77. * O2 majoritaire
  78. LAMBDA = EXP(-1.0D0*TA/T0);
  79. LAMBDA = LAMBDA * (T0 ** (-1.0D0*B));
  80. LAMBDA = LAMBDA * (RYH2I ** (CON-1.0D0)) * RYO2I;
  81. LAMBDA = 2.0D0 * LAMBDA * MH2 * A;
  82. NRYH2I = RYH2I / (1.0D0 + (LAMBDA * DT)) ;
  83. DRYH2 = NRYH2I '-' RYH2I;
  84. DMOL = DRYH2 / MH2 / 2.0D0;
  85. DRYO2 = DMOL * MO2 ;
  86. DRYH2O = -2.0D0 * DMOL * MH2O;
  87. NRYO2I = RYO2I '+' DRYO2 ;
  88. NRYH2OI = RYH2OI '+' DRYH2O ;
  89. 'MESSAGE' ('CHAINE' 'POINT ' &BL1) ;
  90. 'MESSAGE' 'H2 minoritaire ' ;
  91. 'MESSAGE' ('CHAINE' 'DTc = ' (1.0 '/' LAMBDA) ' DT = ' DT) ;
  92. 'FINSI' ;
  93. DEI = (DRYH2 * H0H2) + (DRYO2 * H0O2) + (DRYH2O * H0H2O) ;
  94. DEI = -1.0D0 * DEI ;
  95. 'SINON' ;
  96. NRYO2I = RYO2I ;
  97. NRYH2I = RYH2I ;
  98. NRYH2OI = RYH2OI ;
  99. DEI = 0.0D0 ;
  100. 'FINSI' ;
  101. 'REMPLACER' LIRYH2 &BL1 NRYH2I ;
  102. 'REMPLACER' LIRYO2 &BL1 NRYO2I ;
  103. 'REMPLACER' LIRYH2O &BL1 NRYH2OI ;
  104. 'REMPLACER' LIDE &BL1 DEI ;
  105. 'FIN' BL1 ;
  106.  
  107. NRYH2 = 'MANUEL' 'CHPO' MAIL 1 'H2 ' LIRYH2 ;
  108. NRYO2 = 'MANUEL' 'CHPO' MAIL 1 'O2 ' LIRYO2 ;
  109. NRYH2O = 'MANUEL' 'CHPO' MAIL 1 'H2O ' LIRYH2O ;
  110. DE = 'MANUEL' 'CHPO' MAIL 1 'SCAL' LIDE ;
  111.  
  112.  
  113. 'FINPROC' NRYH2 NRYO2 NRYH2O DE ;
  114.  
  115. ***********************
  116. **** LA TABLE PGAZ ****
  117. ***********************
  118.  
  119. PGAZ = 'TABLE' ;
  120.  
  121. *
  122. **** Ordre des polynoms
  123. *
  124.  
  125. PGAZ . 'NORD' = 4 ;
  126.  
  127. *
  128. **** Especes qui sont dans les equations d'Euler
  129. *
  130.  
  131. PGAZ . 'ESPEULE' = 'MOTS' 'H2 ' 'O2 ' 'H2O ' ;
  132.  
  133. *
  134. **** Espece qui n'y est pas
  135. *
  136.  
  137.  
  138. PGAZ . 'ESPNEULE' = 'N2 ';
  139.  
  140. *
  141.  
  142. PGAZ . 'H2 ' = 'TABLE' ;
  143. PGAZ . 'H2O ' = 'TABLE' ;
  144. PGAZ . 'N2 ' = 'TABLE' ;
  145. PGAZ . 'O2 ' = 'TABLE' ;
  146.  
  147. *
  148. **** R (J/Kg/K)
  149. *
  150.  
  151. PGAZ . 'H2 ' . 'R' = 4130.0 ;
  152. PGAZ . 'H2O ' . 'R' = 461.4 ;
  153. PGAZ . 'N2 ' . 'R' = 296.8 ;
  154. PGAZ . 'O2 ' . 'R' = 259.8 ;
  155.  
  156.  
  157. *
  158. **** Regressions polynomials pour le cvs
  159. *
  160.  
  161. PGAZ . 'H2 ' . 'A' = 'PROG' 9834.91866 0.54273926 0.000862203836
  162. -2.37281455E-07 1.84701105E-11 ;
  163. PGAZ . 'H2O ' . 'A' = 'PROG' 1155.95625 0.768331151 -5.73129958E-05
  164. -1.82753232E-08 2.44485692E-12 ;
  165. PGAZ . 'N2 ' . 'A' = 'PROG' 652.940766 0.288239099 -7.80442298E-05
  166. 8.78233606E-09 -3.05514485E-13 ;
  167. PGAZ . 'O2 ' . 'A' = 'PROG' 575.012333 0.350522002 -0.000128294865
  168. 2.33636971E-08 -1.53304905E-12;
  169.  
  170.  
  171.  
  172. 'SI' FAUX ;
  173.  
  174. *
  175. **** Calcul des enthalpies de formation a 0K
  176. *
  177. * A 298.15K
  178. *
  179. * HOH2 = 0 ;
  180. * HOO2 = 0 ;
  181. * HON2 = 0 ;
  182. * HOH2O = -1.3428D7 J/kg ;
  183. *
  184.  
  185. RH2 = PGAZ . 'H2 ' . 'R' ;
  186. RO2 = PGAZ . 'O2 ' . 'R' ;
  187. RN2 = PGAZ . 'N2 ' . 'R' ;
  188. RH2O = PGAZ . 'H2O ' . 'R' ;
  189.  
  190. A0H2 = 'EXTRAIRE' (PGAZ . 'H2 ' . 'A') 1 ;
  191. A1H2 = 'EXTRAIRE' (PGAZ . 'H2 ' . 'A') 2 ;
  192. A2H2 = 'EXTRAIRE' (PGAZ . 'H2 ' . 'A') 3 ;
  193. A3H2 = 'EXTRAIRE' (PGAZ . 'H2 ' . 'A') 4 ;
  194. A4H2 = 'EXTRAIRE' (PGAZ . 'H2 ' . 'A') 5 ;
  195.  
  196. A0H2O= 'EXTRAIRE' (PGAZ . 'H2O ' . 'A') 1 ;
  197. A1H2O= 'EXTRAIRE' (PGAZ . 'H2O ' . 'A') 2 ;
  198. A2H2O= 'EXTRAIRE' (PGAZ . 'H2O ' . 'A') 3 ;
  199. A3H2O= 'EXTRAIRE' (PGAZ . 'H2O ' . 'A') 4 ;
  200. A4H2O= 'EXTRAIRE' (PGAZ . 'H2O ' . 'A') 5 ;
  201.  
  202. A0N2 = 'EXTRAIRE' (PGAZ . 'N2 ' . 'A') 1 ;
  203. A1N2 = 'EXTRAIRE' (PGAZ . 'N2 ' . 'A') 2 ;
  204. A2N2 = 'EXTRAIRE' (PGAZ . 'N2 ' . 'A') 3 ;
  205. A3N2 = 'EXTRAIRE' (PGAZ . 'N2 ' . 'A') 4 ;
  206. A4N2 = 'EXTRAIRE' (PGAZ . 'N2 ' . 'A') 5 ;
  207.  
  208. A0O2 = 'EXTRAIRE' (PGAZ . 'O2 ' . 'A') 1 ;
  209. A1O2 = 'EXTRAIRE' (PGAZ . 'O2 ' . 'A') 2 ;
  210. A2O2 = 'EXTRAIRE' (PGAZ . 'O2 ' . 'A') 3 ;
  211. A3O2 = 'EXTRAIRE' (PGAZ . 'O2 ' . 'A') 4 ;
  212. A4O2 = 'EXTRAIRE' (PGAZ . 'O2 ' . 'A') 5 ;
  213.  
  214. * On calcule \int_0^298.15 cv(y) dy pour chaque especes
  215.  
  216. T = 298.15 ;
  217. T2 = T * T ;
  218. T3 = T2 * T ;
  219. T4 = T3 * T ;
  220. T5 = T4 * T ;
  221.  
  222. ETHERH2 = ((A0H2 * T) '+' ((A1H2 '/' 2.0) * T2) '+'
  223. ((A2H2 '/' 3.0) * T3) '+' ((A3H2 '/' 4.0) * T4) '+'
  224. ((A4H2 '/' 5.0) * T5)) ;
  225.  
  226. ETHERO2 = ((A0O2 * T) '+' ((A1O2 '/' 2.0) * T2) '+'
  227. ((A2O2 '/' 3.0) * T3) '+' ((A3O2 '/' 4.0) * T4) '+'
  228. ((A4O2 '/' 5.0) * T5)) ;
  229.  
  230. ETHERH2O= ((A0H2O * T) '+' ((A1H2O '/' 2.0) * T2) '+'
  231. ((A2H2O '/' 3.0) * T3) '+' ((A3H2O '/' 4.0) * T4) '+'
  232. ((A4H2O '/' 5.0) * T5)) ;
  233.  
  234. ETHERN2 = ((A0N2 * T) '+' ((A1N2 '/' 2.0) * T2) '+'
  235. ((A2N2 '/' 3.0) * T3) '+' ((A3N2 '/' 4.0) * T4) '+'
  236. ((A4N2 '/' 5.0) * T5)) ;
  237. *
  238.  
  239. ENTHH2 = -1.0D0 * (ETHERH2 '+' (RH2 '*' 298.15)) ;
  240.  
  241. ENTHO2 = -1.0D0 * (ETHERO2 '+' (RO2 '*' 298.15)) ;
  242.  
  243. ENTHN2 = -1.0D0 * (ETHERN2 '+' (RN2 '*' 298.15)) ;
  244.  
  245. ENTHH2O = -1.0D0 * (1.3428D7 '+' ETHERH2O '+' (RH2O '*' 298.15)) ;
  246.  
  247. 'FINSI' ;
  248.  
  249.  
  250. *
  251. **** "Enthalpies" (ou energies) de formations a OK (J/Kg)
  252. *
  253.  
  254. PGAZ . 'H2 ' . 'H0K' = -4.194D6 ;
  255. PGAZ . 'H2O ' . 'H0K' = -1.395D7 ;
  256. PGAZ . 'N2 ' . 'H0K' = -2.953D5 ;
  257. PGAZ . 'O2 ' . 'H0K' = -2.634D5 ;
  258.  
  259.  
  260. ***************************
  261. ***** DOMAINE SPATIAL ****
  262. ***************************
  263.  
  264. *
  265. **** Trois points
  266. *
  267.  
  268. A1 = 0.0D0 0.0D0 ;
  269. DOM1 = 'MANUEL' 'POI1' A1 ;
  270. A2 = 1.0D0 0.0D0 ;
  271. DOM2 = 'MANUEL' 'POI1' A2 ;
  272. A3 = 2.0D0 0.0D0 ;
  273. DOM3 = 'MANUEL' 'POI1' A3 ;
  274.  
  275. DOMTOT = DOM1 'ET' DOM2 'ET' DOM3 ;
  276.  
  277. *
  278. **** FLAM, 'ARRHENIU'
  279. *
  280. * 2 H2 + O2 -> 2 H2O
  281. *
  282. * w = A H(T - Ts) RH2^c RO2 * T^(-b) exp(-Ta/T)
  283. *
  284. * ou
  285. *
  286. * w = taux de reaction MOLAIRE pour la reaction
  287. * chimie, qui est suppose' etre elementaire
  288. *
  289. * H(T - Ts) est la fontion d'Heavyside
  290. *
  291. * N.B. L'energie thermique de chaque epsece est definie par
  292. *
  293. * e = \int_{0}^T cv(y) dy
  294. *
  295.  
  296. Ta = 8310.0;
  297. A = 1.1725D14;
  298. b = 0.91;
  299. con = 2.0;
  300. DT = 0.1;
  301. TS = 'MANUEL' 'CHPO' DOMTOT 1 'SCAL' 1000 'NATURE' 'DISCRET' ;
  302.  
  303. ********************************************
  304. ************* Donnes initiales *************
  305. ********************************************
  306.  
  307. RH2 = PGAZ . 'H2 ' . 'R' ;
  308. RO2 = PGAZ . 'O2 ' . 'R' ;
  309. RN2 = PGAZ . 'N2 ' . 'R' ;
  310. RH2O = PGAZ . 'H2O ' . 'R' ;
  311.  
  312. H0H2 = PGAZ . 'H2 ' . 'H0K' ;
  313. H0O2 = PGAZ . 'O2 ' . 'H0K' ;
  314. H0H2O = PGAZ . 'H2O ' . 'H0K' ;
  315. * H0N2 = PGAZ . 'N2 ' . 'H0K' ;
  316.  
  317. MH2 = 2.016D-3 ;
  318. MO2 = 31.998D-3 ;
  319. MH2O = 18.015D-3 ;
  320. MN2 = 28.013D-3 ;
  321.  
  322. XH2 = 'MANUEL' 'CHPO' DOMTOT 1 'SCAL' ('PROG' 0.3 0.25 0.3) ;
  323.  
  324. XO2 = 'MANUEL' 'CHPO' DOMTOT 1 'SCAL' ('PROG' 0.147 0.147 0.147) ;
  325.  
  326. XH2O = 'MANUEL' 'CHPO' DOMTOT 1 'SCAL' ('PROG' 1.0D-6 1.0D-6 1.0D-6) ;
  327.  
  328. XN2 = ('MANUEL' 'CHPO' DOMTOT 1 'SCAL' 1.0 'NATURE' 'DISCRET' ) '-'
  329. (XH2 '+' XO2 '+' XH2O) ;
  330.  
  331. * Mass des especes dans une mole
  332.  
  333. MH2SMOL = XH2 * MH2 ;
  334. MO2SMOL = XO2 * MO2 ;
  335. MH2OSMOL = XH2O * MH2O ;
  336. MN2SMOL = XN2 * MN2 ;
  337. MTOT = MH2SMOL '+' MO2SMOL '+' MH2OSMOL '+' MN2SMOL ;
  338.  
  339. YH2 = MH2SMOL '/' MTOT ;
  340. YO2 = MO2SMOL '/' MTOT ;
  341. YH2O = MH2OSMOL '/' MTOT ;
  342. YN2 = MN2SMOL '/' MTOT ;
  343.  
  344. RTOT = (YH2 * RH2) '+' (YO2 * RO2) '+' (YH2O * RH2O) '+' (YN2 * RN2) ;
  345.  
  346. * Pression et temperature
  347.  
  348. P = 'MANUEL' 'CHPO' DOMTOT 1 'SCAL' ('PROG' 2.83D6 2.0D6 2.5D6) ;
  349.  
  350. T = 'MANUEL' 'CHPO' DOMTOT 1 'SCAL' ('PROG' 1527D0 2500 1650 ) ;
  351.  
  352. RO = P '/' (RTOT '*' T) ;
  353.  
  354.  
  355. RYH2 = 'NOMC' 'H2 ' (RO * YH2) ;
  356.  
  357. RYO2 = 'NOMC' 'O2 ' (RO * YO2) ;
  358.  
  359. RYH2O = 'NOMC' 'H2O ' (RO * YH2O) ;
  360.  
  361.  
  362. * Dans le cas de detonation dans un tube ou C1 = zone von Neumann et
  363. * C3 = reste du tube, on a
  364. * u ~ 1000 m/s
  365. * gamma ~ 1.25
  366. * donc
  367. * ason ~ 875 m/s -> (u + cson) ~ 1875 m/s
  368. *
  369. * DT ~ DX '/' (u + cson)
  370. *
  371. * DX ~ 0.1 m -> DT ~ 5.00000E-05
  372. *
  373.  
  374. DT = 0.5D-5 ;
  375.  
  376.  
  377. NRYH2 NRYO2 NRYH2O DE = 'FLAM' 'ARRHENIU' TS A b con Ta H0H2 H0O2
  378. H0H2O DT RYH2 RYO2 RYH2O T ;
  379.  
  380. *
  381. **** Conservativite'
  382. *
  383.  
  384. ERRRO = (('EXCO' 'H2 ' NRYH2) '+' ('EXCO' 'O2 'NRYO2) '+'
  385. ('EXCO' 'H2O ' NRYH2O)) '-'
  386. (('EXCO' 'H2 ' RYH2) '+' ('EXCO' 'O2 ' RYO2) '+'
  387. ('EXCO' 'H2O ' RYH2O));
  388.  
  389. ERRRO = ERRRO '/' ('MAXIMUM' RO) ;
  390.  
  391. 'SI' (('MAXIMUM' ERRRO 'ABS') > 1.0D-10) ;
  392. 'ERREUR' 5 ;
  393. 'FINSI' ;
  394.  
  395. NRYH21 NRYO21 NRYH2O1 DE1 = TAUX TS A b con Ta H0H2 H0O2
  396. H0H2O DT RYH2 RYO2 RYH2O T ;
  397.  
  398.  
  399. ERRH2 = ('MAXIMUM' (NRYH2 '-' NRYH21) 'ABS')
  400. '/' ('MAXIMUM' NRYH2 ) ;
  401.  
  402. ERRO2 = 'MAXIMUM' (NRYO2 '-' NRYO21) 'ABS'
  403. '/' ('MAXIMUM' NRYO2 ) ;
  404.  
  405. ERRH2O = 'MAXIMUM' (NRYH2O '-' NRYH2O1) 'ABS'
  406. '/' ('MAXIMUM' NRYH2O ) ;
  407.  
  408.  
  409. ERRDE = 'MAXIMUM' (DE '-' DE1) 'ABS' '/'
  410. ('MAXIMUM' DE) ;
  411.  
  412.  
  413. 'SI' (('MAXIMUM' ('PROG' ERRH2 ERRO2 ERRH2O ERRDE)) > 1.0D-10) ;
  414. 'ERREUR' 5;
  415. 'FINSI' ;
  416.  
  417. *
  418. **** Option HEAVYSIDE: les especes minoritaires brule completement
  419. *
  420.  
  421. NRYH2H NRYO2H NRYH2OH DEH = 'FLAM' 'HEAVYSID' TS H0H2 H0O2
  422. H0H2O RYH2 RYO2 RYH2O T ;
  423.  
  424. *
  425. **** Conservativite'
  426. *
  427.  
  428. ERRRO = (('EXCO' 'H2 ' NRYH2H) '+' ('EXCO' 'O2 'NRYO2H) '+'
  429. ('EXCO' 'H2O ' NRYH2OH)) '-'
  430. (('EXCO' 'H2 ' RYH2) '+' ('EXCO' 'O2 ' RYO2) '+'
  431. ('EXCO' 'H2O ' RYH2O));
  432.  
  433. ERRRO = ERRRO '/' ('MAXIMUM' RO) ;
  434.  
  435. 'SI' (('MAXIMUM' ERRRO 'ABS') > 1.0D-10) ;
  436. 'ERREUR' 5 ;
  437. 'FINSI' ;
  438.  
  439. 'FIN' ;
  440.  
  441.  
  442.  
  443.  
  444.  
  445.  
  446.  

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