C DETO SOURCE CB215821 20/11/25 13:24:54 10792 SUBROUTINE DETO C----------------------------------------------------------------------- C Calcul les conditions Chapman-Jouguet (CJ), Adiabatic Isochore C Complete Combustion (AICC) et Zeldovitch-Neuman-Doringts (ZND) C d'un mélange O2/N2/H2/H2O. C----------------------------------------------------------------------- C C--------------------------- C Phrase d'appel (GIBIANE) : C--------------------------- C C CHP2 CHP3 CHP4 = DETO CHP1 C C------------------------ C Opérandes et résultat : C------------------------ C C CHP1 : CHAMPOINT, composition du mélange, pression et température C composantes : 'O2 ', 'N2 ', 'H2 ', 'H2O ', 'P ', 'T ' C CHP2 : CHAMPOINT, conditions Chapman-Jouguet (CJ) C composantes : 'RCJ ', 'TCJ ', 'PCJ ', 'TAUX', 'VCJ ' C CHP3 : CHAMPOINT, idem Zeldovitch-Neuman-Doringts (ZND) C composantes : 'RZND', 'TZND', 'PZND' C CHP4 : CHAMPOINT, idem Adiabatic Isochore Complete Combustion (AICC) C composantes : 'RAIC', 'TAIC', 'PAIC' C C----------------------------------------------------------------------- C C Langage : ESOPE + FORTRAN77 C C Auteurs : F.DABBENE et E.STUDER 02/99 C C----------------------------------------------------------------------- IMPLICIT INTEGER(I-N) IMPLICIT REAL*8 (A-H,O-Z) CHARACTER*4 LCOMP(6) PARAMETER (NBC=4,NBT=14) DIMENSION AN(NBC),AH2(NBT),AH2O(NBT),AN2(NBT),AO2(NBT) C -INC PPARAM -INC CCOPTIO -INC SMCHPOI -INC SMLMOTS C C- Création du LISTMOTS contenant le nom des composantes du CHAMPOINT : C- 'O2 ', 'N2 ', 'H2 ', 'H2O ', 'P ' et 'T ' décrivant C- respectivement le nombre de moles des différentes especes C- ainsi que la pression et la température du mélange. C IRET = 0 JGN = LOCOMP JGM = 6 SEGINI MLMOTS MOTS(1) = 'H2 ' MOTS(2) = 'O2 ' MOTS(3) = 'H2O ' MOTS(4) = 'N2 ' MOTS(5) = 'P ' MOTS(6) = 'T ' C C----------------------------------------------------------------- C- Evolution des propriétés thermodynamiques des constituants en C- fonction de la température. On retient les modèles de CHEMKIN : C----------------------------------------------------------------- C- cp/R = a1 + a2 T + a3 T**2 + a4 T**3 + a5 T**4 C- H0/RT = a1 + a2 T/2 + a3 T**2 /3 + a4 T**3 /4 + a5 T**4 /5 + a6/T C- S0/R = a1 lnT + a2 T + a3 T**2 /2 + a4 T**3 /3 + a5 T**4 /4 + a7 C C- On a donc besoin de 7 coefficients pour caractériser la dépendance C- à la température. De part et d'autre de la température seuil 1000K C- deux approximations différentes sont proposées -> 14 coefs/espece C C TINT --> Température seuil C IINT --> Nombre de coefficients nécessaire pour T donnée (7) C AH2 --> Coeffs a1 à a14 pour H2 C AH2O --> Coeffs a1 à a14 pour H2O C AN2 --> Coeffs a1 à a14 pour N2 C AO2 --> Coeffs a1 à a14 pour O2 C C Référence : The CHEMKIN Thermodynamic Data Base, R.J.Kee, F.M.Rupley, C J.A.Miller, SANDIA Report, SAND87-8215B. C TINT = 1000.D0 IINT = 7 C AH2(1) = 0.02991423D+02 AH2(2) = 0.07000644D-02 AH2(3) = -0.05633829D-06 AH2(4) = -0.09231578D-10 AH2(5) = 0.01582752D-13 AH2(6) = -0.08350340D+04 AH2(7) = -0.01355110D+02 AH2(8) = 0.03298124D+02 AH2(9) = 0.08249442D-02 AH2(10) = -0.08143015D-05 AH2(11) = -0.09475434D-09 AH2(12) = 0.04134872D-11 AH2(13) = -0.01012521D+05 AH2(14) = -0.03294094D+02 C AH2O(1) = 0.02672146D+02 AH2O(2) = 0.03056293D-01 AH2O(3) = -0.08730260D-05 AH2O(4) = 0.01200996D-08 AH2O(5) = -0.06391618D-13 AH2O(6) = -0.02989921D+06 AH2O(7) = 0.06862817D+02 AH2O(8) = 0.03386842D+02 AH2O(9) = 0.03474982D-01 AH2O(10) = -0.06354696D-04 AH2O(11) = 0.06968581D-07 AH2O(12) = -0.02506588D-10 AH2O(13) = -0.03020811D+06 AH2O(14) = 0.02590233D+02 C AN2(1) = 0.02926640D+02 AN2(2) = 0.01487977D-01 AN2(3) = -0.05684761D-05 AN2(4) = 0.01009704D-08 AN2(5) = -0.06753351D-13 AN2(6) = -0.09227977D+04 AN2(7) = 0.05980528D+02 AN2(8) = 0.03298677D+02 AN2(9) = 0.01408240D-01 AN2(10) = -0.03963222D-04 AN2(11) = 0.05641515D-07 AN2(12) = -0.02444855D-10 AN2(13) = -0.01020900D+05 AN2(14) = 0.03950372D+02 C AO2(1) = 0.03697578D+02 AO2(2) = 0.06135197D-02 AO2(3) = -0.01258842D-05 AO2(4) = 0.01775281D-09 AO2(5) = -0.01136435D-13 AO2(6) = -0.01233930D+05 AO2(7) = 0.03189166D+02 AO2(8) = 0.03212936D+02 AO2(9) = 0.01127486D-01 AO2(10) = -0.05756150D-05 AO2(11) = 0.01313877D-07 AO2(12) = -0.08768554D-11 AO2(13) = -0.01005249D+05 AO2(14) = 0.06034738D+02 C C- Masse molaire en gramme des différents constituants C AN(1) = 2.01594D0 AN(2) = 31.9988D0 AN(3) = 18.01534D0 AN(4) = 28.01340D0 C C- Création du CHPOINT contenant l'ensemble des composantes cherchées, C- le CHAMPOINT étant laissé dans la pile afin que l'opérateur EXTR C- fasse son travail. C CALL ECROBJ('LISTMOTS',MLMOTS) CALL EXCOMP IF (IERR.NE.0) RETURN CALL LIROBJ('CHPOINT',ICHP1,1,IRETOU) IF (IERR.NE.0) RETURN C C- Test du CHAMPOIN : C- 1) le CHPO n'est pas partitionné (i.e. on a toutes les composantes C- en tous les points), sinon erreur C- 110 2 Le %m1:8 de pointeur %i1 n'est pas élémentaire (n<>1) C- 2) les composantes sont dans l'ordre du LISTMOTS, sinon erreur C- 735 2 La composante %m1:4 est inexistante ou mal positionnee C- dans une connectivite C MCHPOI = ICHP1 SEGACT MCHPOI SEGACT MLMOTS N1 = IPCHP(/1) IF (N1.NE.1) THEN MOTERR(1:8) = 'CHAMPOIN' INTERR(1) = ICHP1 CALL ERREUR(110) RETURN ENDIF MSOUPO = IPCHP(1) SEGACT MSOUPO DO 10 I=1,NOCOMP(/2) IF (MOTS(I).NE.NOCOMP(I)) THEN MOTERR= MOTS(I) CALL ERREUR(735) RETURN ENDIF 10 CONTINUE SEGSUP MCHPOI SEGSUP MLMOTS C C- Passage à la phase de calcul C CALL DETO1(MSOUPO,MCHPO1,MCHPO2,MCHPO3, 1 AN,AH2,AH2O,AN2,AO2,TINT,IINT,NBC,NBT) SEGSUP MSOUPO IF (IERR.NE.0) RETURN C C- Ecriture des résultats dans la pile C CALL ECROBJ('CHPOINT',MCHPO3) CALL ECROBJ('CHPOINT',MCHPO2) CALL ECROBJ('CHPOINT',MCHPO1) C RETURN END