Télécharger la plaquette Cast3M


fr en

Télécharger tubedeto2d2.dgibi

Retour à la liste

Numérotation des lignes : 

  1. * fichier :  tubedeto2d2.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4. ****************************************************************
  5. *****      PROPAGATION D UNE DETONATION DANS UN TUBE       *****
  6. *****      MODELE COMBUSTION H2-air de PLEXUS              *****
  7. *****      A. BECCANTINI, SFME/LTMF, 15.01.02              *****
  8. ****************************************************************
  9. *
  10. * The file is divided into 5 parts
  11. *
  12. * 1) mesh
  13. * 2) initial conditions and gas properties
  14. * 3) parameters used in the computation
  15. * 4) the main part (where the computation is realised)
  16. * 5) the post-treatment
  17. *
  18.  
  19.    'OPTION' 'ECHO' 1 'DIME' 2
  20.          'TRAC' 'X' ;
  21.  
  22. ****************************************************************
  23. ****************************************************************
  24. *****      1D mesh                                         *****
  25. ****************************************************************
  26. ****************************************************************
  27. *
  28. *
  29. *  A6 ---- A5 ------------------- A4 ----------------------- A8
  30. *   |      |                       |                         |
  31. *   |  AR  |        ZONE1          |         ZONE2           |
  32. *   |      |                       |                         |
  33. *  A1 ---- A2 ------------------- A3 ----------------------- A7
  34. *   |      |                       |                         |
  35. *   |  L1  |          L2           |           L3            |
  36. *   |<---->|<--------------------->|<----------------------->|
  37. *
  38. *  AR = activation region
  39. *  ZONE1
  40. *  ZONE2
  41. *
  42.  
  43.  RAF = 1  ;
  44.  
  45.  NL2SL1 = 49 ;
  46.  NL3SL1 = 100 ;
  47.  L1 = 0.4 ;
  48.  
  49.  L2 = NL2SL1 '*' L1 ;
  50.  L3 = NL3SL1 '*' L1 ;
  51.  
  52. * If  NL2SL1 = 49 and NL3SL1 = 50, then L2 = 19.6, L3 = 20.0 ;
  53.  
  54.  N1 = RAF ;
  55.  N2 = NL2SL1 '*' N1 ;
  56.  N3 = NL3SL1 '*' N1 ;
  57.  
  58. * For sake of simplicity, we will only consider one mesh in y-direction
  59.  
  60.  DX = L1 '/' N1 ;
  61.  NY = 4 ;
  62.  NZ = NY ;
  63.  H1 = NY '*' DX ;
  64.  
  65.  A1 = 0.0 0.0 ;
  66.  A2 = L1  0.0 ;
  67.  A3 = (L1 '+' L2) 0.0 ;
  68.  A4 = (L1 '+' L2) H1 ;
  69.  A5 = L1  H1 ;
  70.  A6 = 0.0 H1 ;
  71.  A7 = (L1 '+' L2 '+' L3) 0.0 ;
  72.  A8 = (L1 '+' L2 '+' L3) H1 ;
  73.  
  74.  'OPTION' 'ELEM' 'SEG2' ;
  75.  A1A2 = A1 'DROIT' N1 A2 ;
  76.  A2A3 = A2 'DROIT' N2 A3 ;
  77.  A3A4 = A3 'DROIT' NY A4 ;
  78.  A4A5 = A4 'DROIT' N2 A5 ;
  79.  A5A2 = A5 'DROIT' NY A2 ;
  80.  A2A5 = 'INVERSE' A5A2 ;
  81.  A5A6 = A5 'DROIT' N1 A6 ;
  82.  A6A1 = A6 'DROIT' NY A1 ;
  83.  A3A7 = A3 'DROIT' N3 A7 ;
  84.  A7A8 = A7 'DROIT' NY A8 ;
  85.  A8A4 = A8 'DROIT' N3 A4 ;
  86.  A4A3 = 'INVERSE' A3A4 ;
  87.  
  88. *
  89. **** 2D mesh
  90. *
  91.  
  92. * DOM1 is the activation region
  93. * DOM2 is zone 1
  94. * DOM3 is zone 2
  95.  
  96.  'OPTION' 'ELEM' 'QUA4' ;
  97.  DOM1 = 'TRANSLATION'  A2A5 N1 ((-1.0 '*' L1) 0.0) ;
  98.  DOM2 = 'TRANSLATION'  A2A5 N2 (L2 0.0) ;
  99.  DOM3 = 'TRANSLATION'  A3A4 N3 (L3 0.0) ;
  100.  
  101. * DOMTOT = total region
  102.  
  103.  DOMTOT = DOM1 'ET' DOM2 'ET' DOM3 ;
  104.  'ELIMINATION' (1.0D-3 '/' RAF) DOMTOT ;
  105.  
  106. *
  107. **** The tables 'DOMAINE'
  108. *
  109.  
  110.  $DOMTOT = 'MODELISER' DOMTOT 'EULER' ;
  111.  $DOM1   = 'MODELISER' DOM1 'EULER' ;
  112.  $DOM2   = 'MODELISER' DOM2 'EULER' ;
  113.  $DOM3   = 'MODELISER' DOM3 'EULER' ;
  114.  
  115.  TDOMTOT = 'DOMA' $DOMTOT 'VF' ;
  116.  TDOM1   = 'DOMA' $DOM1   'VF' ;
  117.  TDOM2   = 'DOMA' $DOM2   'VF' ;
  118.  TDOM3   = 'DOMA' $DOM3   'VF' ;
  119.  
  120.  QDOMTOT = TDOMTOT  . 'QUAF' ;
  121.  QDOM1   = TDOM1    . 'QUAF' ;
  122.  QDOM2   = TDOM2    . 'QUAF' ;
  123.  QDOM3   = TDOM3    . 'QUAF' ;
  124.  
  125.  'ELIMINATION' QDOMTOT (0.001 '/' RAF) QDOM1 ;
  126.  'ELIMINATION' QDOMTOT (0.001 '/' RAF) QDOM2 ;
  127.  'ELIMINATION' QDOMTOT (0.001 '/' RAF) QDOM3 ;
  128.  
  129. *
  130. **** MOD1 = model (created just to display the CHAMELEMs)
  131. *
  132.  
  133.  MOD1 =  'MODELISER' ('DOMA' $DOMTOT 'MAILLAGE') 'THERMIQUE' ;
  134.  
  135. *
  136. **** Line for graphics
  137. *
  138.  
  139.  'OPTION' 'ELEM' 'SEG2' ;
  140.  P1 = (DX '/' 2) (DX '/' 2);
  141.  P2 = ((L1 '+' L2 '+' L3) '-' (DX '/' 2)) (DX '/' 2);
  142.  LIGEVO = (P1 'DROIT' (N1 '+' N2 '+' N3 '-' 1) P2) 'COULEUR' 'ROUG'  ;
  143.  'ELIMINATION' LIGEVO ('DOMA' $DOMTOT 'CENTRE') (0.001 '/' RAF) ;
  144.  
  145. * 'OPTION' 'SAUVER' ('CHAINE' 'mesh.sauv_raf' RAF) ;
  146. * 'SAUVER' RAF DOMTOT $DOMTOT $DOM1 $DOM2 $DOM3 LIGEVO MOD1 ;
  147.  
  148. ****************************************************************
  149. ****************************************************************
  150. *****      Initial conditions and gas properties.          *****
  151. ****************************************************************
  152. ****************************************************************
  153.  
  154. *
  155. **** 0 from a numerical point of view
  156. *
  157.  
  158.  zero = 1.0d-9 ;
  159.  
  160. *
  161. **** (Non-Uniform) initial conditions.(DOM1 et DOM2)
  162. *    We compute the conservative variables.
  163. *
  164. *
  165. * Contact discontinuity between DOM2 and DOM3 (same pressure,
  166. * different temperature)
  167. *
  168.  
  169.  PINI = ('MANUEL' 'CHPO' (('DOMA' $DOM1 'CENTRE') 'ET'
  170.         ('DOMA' $DOM2 'CENTRE')) 1 'SCAL' 99700 'NATU' 'DISCRET') 'ET'
  171.         ('MANUEL' 'CHPO' (('DOMA' $DOM3 'CENTRE')) 1 'SCAL' 99700
  172.         'NATU' 'DISCRET') ; 
  173.  TINI = ('MANUEL' 'CHPO' (('DOMA' $DOM1 'CENTRE') 'ET'
  174.         ('DOMA' $DOM2 'CENTRE')) 1 'SCAL' 285. 'NATU' 'DISCRET') 'ET'
  175.         ('MANUEL' 'CHPO' (('DOMA' $DOM3 'CENTRE')) 1 'SCAL' 400
  176.         'NATU' 'DISCRET') ; 
  177.  WINI = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 2 'UX' 0.0 'UY' 0.0;
  178.  
  179. * We suppose that the combustion is complete everywhere, i.e. the
  180. * progress variable is 1 after the detonation transition.
  181.  
  182. * CSIMAX = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'SCAL' 1.0 ;
  183.  
  184. *
  185. **** From initial conditions on molar fractions, we duduce
  186. *    the data for the virtual mixture of burned/unburned
  187. *    gases
  188.  
  189. *
  190. **** Molar fractions (before combustion)
  191. *
  192.  
  193.  XH21 = ('MANUEL' 'CHPO' (('DOMA' $DOM1 'CENTRE') 'ET'
  194.         ('DOMA' $DOM2 'CENTRE')) 1 'H2' 0.3 'NATU' 'DISCRET') 'ET'
  195.         ('MANUEL' 'CHPO' (('DOMA' $DOM3 'CENTRE')) 1 'H2' 0.2
  196.          'NATU' 'DISCRET') ;
  197.  XH2O1 = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'H2O' zero
  198.          'NATURE' 'DISCRET' ;
  199.  XO21 = ('MANUEL' 'CHPO' (('DOMA' $DOM1 'CENTRE') 'ET'
  200.         ('DOMA' $DOM2 'CENTRE')) 1 'O2' 0.147 'NATU' 'DISCRET') 'ET'
  201.         ('MANUEL' 'CHPO' (('DOMA' $DOM3 'CENTRE')) 1 'O2' 0.2
  202.          'NATU' 'DISCRET') ;
  203.  XN = XH21 'ET' XH2O1 'ET' XO21 ;
  204.  CHUN = ('MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'SCAL' 1.0D0) ;
  205.  XN21 =  CHUN '-' ('PSCAL' XN CHUN  ('MOTS' 'H2' 'O2' 'H2O')
  206.          ('MOTS' 'SCAL' 'SCAL' 'SCAL') ) ;
  207.  XN21 = 'NOMC' 'N2' XN21 'NATURE' 'DISCRET' ;
  208.  XN = XN 'ET' XN21 ;
  209.  
  210. *
  211. ***** Gas properties
  212. *
  213.  
  214.  PGAZ = 'TABLE' ;
  215.  
  216. * Polynomial degree of specific heats
  217.  
  218.  PGAZ . 'NORD' = 4 ;
  219.  
  220. * Species explicitly treated in the Euler Equations
  221.  
  222.  PGAZ . 'ESPEULE' = 'MOTS' 'H2  ' 'O2  ' 'H2O ' ;
  223.  
  224. * Species non explicitly treated
  225.  
  226.  PGAZ . 'ESPNEULE' = 'N2  ';
  227.  
  228. * Single gases properties
  229.  
  230.  PGAZ .  'H2  ' = 'TABLE'  ;
  231.  PGAZ .  'H2O ' = 'TABLE'  ;
  232.  PGAZ .  'N2  ' = 'TABLE'  ;
  233.  PGAZ .  'O2  ' = 'TABLE'  ;
  234.  
  235. * R (J/Kg/K)
  236.  
  237.  mH2 = 2. '*' 1.00797E-3 ;
  238.  mo2 = 2. '*' 15.9994E-3 ;
  239.  mH2O = mh2 '+' (0.5 '*' mo2) ;
  240.  mN2 = 2 '*' 14.0067E-3 ;
  241.  RGAS = 8.31441 ;
  242.  
  243.  PGAZ .  'H2  ' . 'R' = RGAS '/' mh2 ;
  244.  PGAZ .  'H2O ' . 'R' = RGAS '/' mh2o ;
  245.  PGAZ .  'N2  ' . 'R' = RGAS '/' mn2 ;
  246.  PGAZ .  'O2  ' . 'R' = RGAS '/' mo2 ;
  247.  
  248. * Polynomials regressions coefficients
  249.  
  250.  PGAZ .  'H2  ' . 'A' = 'PROG'  9834.91866 0.54273926 0.000862203836
  251.                                -2.37281455E-07 1.84701105E-11 ;
  252.  PGAZ .  'H2O ' . 'A' = 'PROG' 1155.95625 0.768331151 -5.73129958E-05
  253.                               -1.82753232E-08 2.44485692E-12 ;
  254.  PGAZ .  'N2  ' . 'A' = 'PROG' 652.940766 0.288239099 -7.80442298E-05
  255.                               8.78233606E-09 -3.05514485E-13 ;
  256.  PGAZ .  'O2  ' . 'A' = 'PROG' 575.012333  0.350522002 -0.000128294865
  257.                               2.33636971E-08 -1.53304905E-12;
  258.  
  259. * Formation enthalpies at 0K (J/Kg))
  260. *
  261. *     h_i(0K)  = h_i(T0) '-' \int_0^{T0} cp_i(x) dx 
  262. *              = h_i(T0) '-' (\int_0^{T0} cv_i(x) dx '+' R_i * T0)
  263. *     
  264.  
  265.  PGAZ .  'H2  ' . 'H0K' = -4.195D6 ;
  266.  PGAZ .  'H2O ' . 'H0K' = -1.395D7 ;
  267.  PGAZ .  'N2  ' . 'H0K' = -2.953D5 ;
  268.  PGAZ .  'O2  ' . 'H0K' = -2.634D5 ;
  269.  
  270.  
  271. * 'OPTION' 'SAUVER' 'gas_properties.sauv' ;
  272. * 'SAUVER' PGAZ ;
  273.  
  274. * 
  275. **** We compute the mass fractions
  276. *
  277.  
  278.  MMOL = ('MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'H2' Mh2
  279.          'NATURE' 'DISCRET') 'ET'
  280.         ('MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'O2' Mo2
  281.          'NATURE' 'DISCRET') 'ET'
  282.         ('MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'H2O' Mh2o
  283.          'NATURE' 'DISCRET') 'ET'
  284.         ('MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'N2' Mn2
  285.          'NATURE' 'DISCRET') ;
  286.  
  287.  MTOT = 'PSCAL' MMOL XN ('MOTS' 'H2' 'O2' 'H2O' 'N2')
  288.    ('MOTS' 'H2' 'O2' 'H2O' 'N2') ;
  289.  
  290.  YH21 = (XH21 '*' MH2) '/' MTOT ;
  291.  YO21 = (XO21 '*' MO2) '/' MTOT ;
  292.  YN21 = (XN21 '*' MN2) '/' MTOT ;
  293.  YH2O1 = (XH2O1 '*' MH2O) '/' MTOT ;
  294.  
  295. *
  296. **** The conservative variables
  297. *
  298.  
  299.  RINI = (('NOMC' 'SCAL' yh21) '*' (PGAZ .'H2' . 'R')) '+'
  300.         (('NOMC' 'SCAL' yo21) '*' (PGAZ .'O2' . 'R')) '+'
  301.         (('NOMC' 'SCAL' yh2o1) '*' (PGAZ .'H2O' . 'R')) '+'
  302.         (('NOMC' 'SCAL' yn21) '*' (PGAZ .'N2' . 'R'))
  303.    ;
  304.  
  305.  RN = PINI '/' (RINI '*' TINI) ;
  306.  
  307. * Internal energy
  308.  
  309.  eini = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'SCAL' 0.0 ;
  310.  Tfois = Tini ;
  311.  
  312.  'REPETER' BL1 ((PGAZ . 'NORD') '+' 1) ;
  313.  
  314.     acel = (('NOMC' 'SCAL' yh21) '*' ('EXTRAIRE' (PGAZ .'H2  ' . 'A')
  315.             &BL1)) '+'
  316.            (('NOMC' 'SCAL' yo21) '*' ('EXTRAIRE' (PGAZ .'O2  ' . 'A')
  317.             &BL1)) '+'
  318.            (('NOMC' 'SCAL' yh2o1) '*' ('EXTRAIRE' (PGAZ .'H2O ' . 'A')
  319.             &BL1)) '+'
  320.            (('NOMC' 'SCAL' yn21) '*' ('EXTRAIRE' (PGAZ .'N2  ' . 'A')
  321.             &BL1))  ;
  322.  
  323.     eini = eini '+' (acel '*' (Tfois '/' &BL1)) ;
  324.  
  325.     Tfois = Tfois * tini ;
  326.  
  327.  'FIN' BL1 ;
  328.  
  329.  RETN = RN '*' (eini '+' (0.5 '*' ('PSCAL' WINI WINI
  330.    ('MOTS' 'UX' 'UY') ('MOTS' 'UX' 'UY')))) ;
  331.  GN = RN '*' WINI ;
  332.  YN = YH21 '+' YO21 '+' YH2O1 ;
  333.  RYN = RN '*' YN ;
  334.  
  335.  'SI' VRAI ;
  336. *   Test  
  337.     VN P T Y GAMN = 'PRIM' 'PERFTEMP' PGAZ RN GN RETN RYN ;
  338.     PCEL = 'MAXIMUM' PINI ;
  339.     TCEL = 'MAXIMUM' TINI ;
  340.     'SI' (('MAXIMUM' ((P '-' pini) '/' PCEL) 'ABS')> 1.0D-3) ;
  341.        'MESSAGE' 'Probleme calcul energie totale' ;
  342.        'ERREUR' 21 ;
  343.     'FINSI' ;
  344.     'SI' (('MAXIMUM' ((T '-' Tini) '/' TCEL)) > 1.0D-3) ;
  345.        'MESSAGE' 'Probleme calcul energie totale' ;
  346.        'ERREUR' 21 ;
  347.     'FINSI' ;
  348.  'FINSI' ;
  349.  
  350. ****************************************************************
  351. ****************************************************************
  352. *****           Parameters for the computations           ******
  353. ****************************************************************
  354. ****************************************************************
  355.  
  356. * Upwind scheme
  357. * METO = 'VLH' ;
  358.  METO = 'SS' ;
  359.  
  360. * Iterations
  361. * Final time
  362. * Safety Factor fot the time step
  363. * Second order reconstruction?
  364.  
  365.  NITER = 10000 ;
  366.  TFINAL = 5.0D-3 ;
  367.  SAFFAC = 0.7 ;
  368.  LOGSO = VRAI ;
  369.  
  370. * 'OPTION' 'SAUVER' 'parameters.sauv' ;
  371. * 'SAUVER'  METO  NITER  TFINAL  SAFFAC  LOGSO ;
  372.  
  373. ****************************************************************
  374. ****************************************************************
  375. *****           The computation                           ******
  376. ****************************************************************
  377. ****************************************************************
  378.  
  379. ************************
  380. **** PROCEDURE DETO ****
  381. ************************
  382.  
  383.  'DEBPROC'  DETO ;
  384.  'ARGUMENT' PGAZ*'TABLE' TN*'CHPOINT' RYN*'CHPOINT' DTI*'FLOTTANT' ;
  385.  
  386. *
  387. * dRYO2/dT  = MO2 * A * (T ** -B) EXP( - TA / T)
  388. *             (RYH2 ** C) * (RYO2)
  389. *             
  390. **** Les parametres physiques
  391. *
  392. *    Ts = temperature de seuil
  393. *
  394.  
  395.     TS = 'MANUEL' 'CHPO' ('EXTRAIRE' TN 'MAILLAGE') 1 'SCAL' 600 ;
  396.     TA = 8310.;
  397.     A = 1.1725D14;
  398.     B = 0.91;
  399.     CON = 2;
  400.  
  401.     H0H2  = PGAZ .  'H2  ' . 'H0K' ;
  402.     H0O2  = PGAZ .  'O2  ' . 'H0K' ;
  403.     H0H2O = PGAZ .  'H2O ' . 'H0K' ;
  404.  
  405.     RYH2 = 'EXCO' 'H2' RYN 'H2' 'NATU' 'DISCRET' ;
  406.  
  407.     RYO2 = 'EXCO' 'O2' RYN 'O2' 'NATU' 'DISCRET' ;
  408.  
  409.     RYH2O = 'EXCO' 'H2O' RYN 'H2O' 'NATU' 'DISCRET' ;
  410.  
  411.     RYH2 RYO2 RYH2O DELTAE = 'FLAM' 'ARRHENIU' TS A B CON TA
  412.                             H0H2 H0O2 H0H2O DTI RYH2 RYO2 RYH2O
  413.                             TN ;
  414.  
  415.     RY =  RYH2 'ET' RYO2 'ET' RYH2O ;
  416.  
  417.  'FINPROC' RY DELTAE ;
  418.  
  419. ****************************
  420. **** FIN PROCEDURE DETO ****
  421. ****************************
  422.  
  423. *
  424. **** Chemical time step 
  425. *
  426. *    We compute d(RHY2)/dT at the von Neumann state
  427. *
  428. *    Molar reaction rate
  429. *
  430.  
  431.  TVN = 2000 ;
  432.  TA = 8310.;
  433.  A = 1.1725D14;
  434.  B = 0.91;
  435.  CON = 2;
  436.  
  437.  H0H2  = PGAZ .  'H2  ' . 'H0K' ;
  438.  H0O2  = PGAZ .  'O2  ' . 'H0K' ;
  439.  H0H2O = PGAZ .  'H2O ' . 'H0K' ;
  440.  
  441.  RYNRED = 'REDU' RYN ('DOMA' $DOM2 'CENTRE') ; 
  442.  RYH2 = 'MAXIMUM' ('EXCO' 'H2' RYNRED 'SCAL') ;
  443.  RYO2 = 'MAXIMUM' ('EXCO' 'O2' RYNRED 'SCAL') ;
  444.  
  445.  OMEGA = A '*' (RYH2 '**' CON) '*' RYO2 '*' (TVN '**' (-1.0 '*' B))
  446.    '*' ('EXP' (-1.0 '*' TA '*' (TVN '**' -1))) ;
  447.  
  448. * OMEGA in (mol/m3/s)
  449. * OMEGAH2 in (kg/m3/s)
  450.  
  451.  OMEGAH2 =  OMEGA '*' 2 '*' 2E-3 ;
  452.  DELTATC = RYH2 '/' OMEGAH2 ; 
  453.  
  454. * NB DELTATC << DELTAT_CON
  455.  
  456. *
  457. **** Zone d activation = etat AICC
  458. *
  459.  
  460.  VN PN TN YN GAMN = 'PRIM' 'PERFTEMP' PGAZ RN GN RETN RYN ;
  461.  RYH2 = ('EXCO' 'H2' RYN 'H2' 'NATU' 'DISCRET') ;
  462.  RYO2 = ('EXCO' 'O2' RYN 'O2' 'NATU' 'DISCRET') ;
  463.  RYH2O = ('EXCO' 'H2O' RYN 'H2O' 'NATU' 'DISCRET') ;
  464.  TN = TN '-' ('REDU' TN ('DOMA' $DOM1 'CENTRE')) '+'
  465.       ('MANUEL' 'CHPO' ('DOMA' $DOM1 'CENTRE') 1 'SCAL' TVN) ;
  466.  TS = 'MANUEL' 'CHPO' ('EXTRAIRE' TN 'MAILLAGE') 1 'SCAL' (0.9 '*' TVN)
  467.     ;
  468.  RYH2 RYO2 RYH2O DELTARE = 'FLAM' 'ARRHENIU' TS A B CON TA
  469.                             H0H2 H0O2 H0H2O (1.0D5 '*' DELTATC) 
  470.                             RYH2 RYO2 RYH2O TN ;
  471.  
  472.  RYN =  RYH2 'ET' RYO2 'ET' RYH2O ;
  473.  RETN = RETN '+' DELTARE ;
  474.  
  475.  RN0 = 'COPIER' RN ;
  476.  GN0 = 'COPIER' GN ;
  477.  RETN0 = 'COPIER' RETN ;
  478.  RYN0 = 'COPIER' RYN ;
  479.  
  480. *
  481. **** Parameter for the time loop 
  482. *
  483.  
  484. * Names of the residuum components
  485.  
  486.  LISTINCO = 'MOTS' 'RN' 'RUX' 'RUY' 'RETN' 'H2' 'O2' 'H2O' ;
  487.  
  488. *
  489. **** Geometrical coefficient to compute gradients
  490. *
  491.  
  492.  GRADR CACCA COEFSCAL =  'PENT' $DOMTOT 'CENTRE' 'EULESCAL' 'NOLIMITE'
  493.            ('MOTS' 'SCAL') RN ;
  494.  
  495.  GRADV CACCA COEFVECT =  'PENT' $DOMTOT 'CENTRE' 'EULEVECT' 'NOLIMITE'
  496.           ('MOTS' 'UX' 'UY') GN ;
  497.  
  498.  TPS = 0.0 ;
  499.  
  500. *
  501. **** The chemistry time step is too little
  502. *    We take DT_CHEM > DELTATC 
  503. *
  504. *
  505.  
  506.  DT_CHEM = 50 '*' DELTATC ;
  507.  
  508. *
  509. **** Temporal loop
  510. *
  511.  
  512.  'MESSAGE' ;
  513.  'MESSAGE' ('CHAINE' 'Methode = ' METO) ;
  514.  'MESSAGE' ('CHAINE' 'SAFFAC      = ' SAFFAC) ;
  515.  'MESSAGE' ;
  516.  
  517.   VN PN TN YN SN GAMMAN = 'PRIM' 'PERFTEMP' PGAZ RN GN RETN RYN RSN ;
  518.  
  519.   TNM1 = 'COPIER' TN ;
  520.  
  521.  'TEMPS' 'ZERO' ;
  522.  'REPETER' BL1 NITER ;
  523.  
  524. *
  525. **** Primitive variables
  526. *
  527.  
  528.     VN PN TN YN GAMMAN = 'PRIM' 'PERFTEMP' PGAZ RN GN RETN RYN 
  529.       TNM1  ;
  530.  
  531.     TNM1 = 'COPIER' TN ;
  532.  
  533.     'SI' LOGSO ;
  534.  
  535.        GRADR ALR = 'PENT' $DOMTOT 'CENTRE' 'EULESCAL' 'LIMITEUR'
  536.           ('MOTS' 'SCAL') RN  'GRADGEO' COEFSCAL ;
  537.  
  538.        GRADP ALP = 'PENT' $DOMTOT 'CENTRE' 'EULESCAL' 'LIMITEUR'
  539.            ('MOTS' 'SCAL') PN  'GRADGEO' COEFSCAL ;
  540.  
  541.        GRADV ALV = 'PENT' $DOMTOT 'CENTRE' 'EULEVECT' 'LIMITEUR'
  542.           ('MOTS' 'UX' 'UY') VN   'GRADGEO'  COEFVECT ;
  543.  
  544.        GRADY ALY = 'PENT' $DOMTOT 'CENTRE' 'EULESCAL' 'LIMITEUR'
  545.            ('MOTS' 'H2' 'O2' 'H2O') YN  'GRADGEO' COEFSCAL ;
  546.  
  547.        ROF VITF PF YF = 'PRET' 'PERFTEMP'  2 1
  548.                             $DOMTOT  PGAZ
  549.                             RN    GRADR  ALR
  550.                             VN    GRADV  ALV
  551.                             PN    GRADP  ALP
  552.                             YN    GRADY  ALY ;
  553.  
  554.     'SINON' ;
  555.  
  556.        ROF VITF PF YF  = 'PRET' 'PERFTEMP'  1 1 
  557.                               $DOMTOT  PGAZ
  558.                               RN   
  559.                               VN   
  560.                               PN   
  561.                               YN    ;
  562.  
  563.     'FINSI' ;
  564.  
  565.     RESIDU DELTAT =  'KONV' 'VF' 'PERFTEMP' 'RESI' METO
  566.          $DOMTOT PGAZ LISTINCO  ROF VITF PF YF ;
  567.  
  568.     DT_CON = SAFFAC '*' DELTAT ;     
  569. *
  570. **** The time step linked to tps
  571. *
  572.  
  573.     DTTPS = (TFINAL '-' TPS) * (1. '+' ZERO) ;
  574.  
  575. *
  576. **** Total time step
  577. *
  578.  
  579.     DTMIN = 'MINIMUM' ('PROG' DT_CON DTTPS DT_CHEM) ;
  580.  
  581. *
  582. **** Increment of the variables (convection)
  583. *
  584.  
  585.     RESIDU = DTMIN '*' RESIDU ;
  586.  
  587.     DRN = 'EXCO' 'RN' RESIDU 'SCAL' ;
  588.     DGN = 'EXCO' ('MOTS' 'RUX' 'RUY') RESIDU ('MOTS' 'UX' 'UY') ;
  589.     DRETN = 'EXCO' 'RETN' RESIDU 'SCAL' ;
  590.     DRYN  = 'EXCO' ('MOTS' 'H2' 'O2' 'H2O') RESIDU
  591.             ('MOTS' 'H2' 'O2' 'H2O') ;
  592.  
  593.  
  594.     TPS  = TPS '+' DTMIN ;
  595.     RN   = RN '+' DRN ;
  596.     GN   = GN '+' DGN ;
  597.     RETN = RETN '+' DRETN ;
  598.     RYN = RYN '+' DRYN ;
  599.  
  600. *
  601. **** Increment of the variables (chemical source)
  602. *
  603.  
  604.     VN PN TN YN GAMMAN = 'PRIM' 'PERFTEMP' PGAZ RN GN RETN RYN 
  605.       TNM1  ;
  606.  
  607.     RYN DRETN = DETO PGAZ TN RYN DTMIN ;
  608.     RETN = RETN '+' DRETN ;
  609.  
  610.  
  611.     'SI' (((&BL1 '/' 20) '*' 20) 'EGA' &BL1) ;
  612.         'MESSAGE' ('CHAINE' 'ITER =' &BL1 '  TPS =' TPS) ;
  613.     'FINSI' ;
  614.  
  615.     'SI' (TPS > TFINAL) ;
  616.        'QUITTER' BL1 ;
  617.     'FINSI' ;
  618.  
  619.  'FIN' BL1 ;
  620.  'TEMPS' ;
  621.  
  622. * 'OPTION' 'SAUVER' ('CHAINE' 'result.sauv_' RAF 'tps_' TPS ) ;
  623. * 'SAUVER'  ;
  624.  
  625. ****************************************************************
  626. ****************************************************************
  627. *****           The post-treatment                        ******
  628. ****************************************************************
  629. ****************************************************************
  630.  
  631. * 'OPTI' 'REST' 'result.sauv_1tps_5' ;
  632. * 'REST' ;
  633.  
  634.  GRAPH = VRAI ;
  635.  GRAPH = FAUX ;
  636.  
  637. * 
  638. **** The mesh
  639. *
  640.  
  641.  'SI' GRAPH ;
  642.     'TRACER' DOMTOT 'TITR' ('CHAINE' 'Maillage: '
  643.        ('NBEL' DOMTOT) ' elts');
  644.  'FINSI' ;
  645.  
  646. *
  647. **** Initial conditions
  648. *
  649.  
  650.   MOD1 =  'MODELISER'  ('DOMA' $DOMTOT 'MAILLAGE') 'THERMIQUE' ;
  651.  
  652.  'SI' GRAPH ;
  653.     VN PN TN YN GAMMAN = 'PRIM' 'PERFTEMP'
  654.        PGAZ RN0 GN0 RETN0 RYN0  ;
  655.     CHM_RN  =  'KCHA'  $DOMTOT  'CHAM'  RN0  ;
  656.     CHM_VN  =  'KCHA'  $DOMTOT  'CHAM'  VN  ;
  657.     CHM_PN  =  'KCHA'  $DOMTOT  'CHAM'  PN  ;
  658.     CHM_TN  =  'KCHA'  $DOMTOT  'CHAM'  TN  ;
  659.     CHM_YN  =  'KCHA'  $DOMTOT  'CHAM'  YN  ;
  660.     'TRAC' CHM_RN  MOD1 
  661.        'TITR' ('CHAINE' 'rho at t=' 0.0)   ;
  662.     'TRAC' CHM_VN  MOD1 
  663.        'TITR' ('CHAINE' 'v at t= '  0.0)   ;
  664.     'TRAC' CHM_PN  MOD1 
  665.        'TITR' ('CHAINE' 'p at t= '  0.0)   ;
  666.     'TRAC' CHM_TN  MOD1 
  667.        'TITR' ('CHAINE' 'T at t= '  0.0)   ;
  668.     'TRAC' CHM_YN  MOD1 
  669.        'TITR' ('CHAINE' 'y at t= '  0.0) ;
  670.  'FINSI' ;
  671.  
  672. *
  673. **** Results
  674. *
  675.  
  676.   VN PN TN YN  GAMMAN = 'PRIM' 'PERFTEMP' PGAZ RN GN RETN RYN ;
  677.  
  678.  'SI' GRAPH ;
  679.  
  680.     CHM_RN  =  'KCHA'  $DOMTOT  'CHAM'  RN  ;
  681.     CHM_VN  =  'KCHA'  $DOMTOT  'CHAM'  VN  ;
  682.     CHM_PN  =  'KCHA'  $DOMTOT  'CHAM'  PN  ;
  683.     CHM_TN  =  'KCHA'  $DOMTOT  'CHAM'  TN  ;
  684.     CHM_YN  =  'KCHA'  $DOMTOT  'CHAM'  YN  ;
  685.  
  686.     'TRAC' CHM_RN  MOD1 
  687.        'TITR' ('CHAINE' 'rho at t=' TPS)   ;
  688.     'TRAC' CHM_VN  MOD1 
  689.        'TITR' ('CHAINE' 'v at t= '  TPS)   ;
  690.     'TRAC' CHM_PN  MOD1 
  691.        'TITR' ('CHAINE' 'p at t= '  TPS)   ;
  692.     'TRAC' CHM_TN  MOD1 
  693.        'TITR' ('CHAINE' 'T at t= '  TPS)   ;
  694.     'TRAC' CHM_YN  MOD1 
  695.        'TITR' ('CHAINE' 'y at t= '  TPS)   ;
  696.  
  697.  'FINSI' ;
  698.  
  699. *
  700. *** Evolution Objects
  701. *
  702.  
  703.  SS = 'COORDONNEE' 1 LIGEVO;
  704.  
  705.  evxx = 'EVOL' 'CHPO' ss LIGEVO ;
  706.  lxx = 'EXTRAIRE'  evxx 'ORDO' ;
  707.  
  708.  x0 = 'MINIMUM' lxx;
  709.  x1 = 'MAXIMUM' lxx ;
  710.  
  711.  titolo  = 'CHAINE'  ' at t = ' TPS '  Methode = ' METO
  712.    ' SF = ' SAFFAC ;
  713.  
  714.  VNN = 'EXCO' 'UX' VN ;
  715.  VNT = 'EXCO' 'UY' VN ;
  716.  
  717. * rho
  718.  
  719.  evro = 'EVOL'  'CHPO' RN LIGEVO ;
  720.  lro  = 'EXTRAIRE'  evro 'ORDO' ;
  721.  evro = 'EVOL' 'MANU' 'x' lxx 'ro' lro ;
  722.  tro  = 'CHAINE' 'ro ' titolo ;
  723.  
  724. * u
  725.  
  726.  evu  = 'EVOL' 'CHPO' VNN LIGEVO ;
  727.  lu   = 'EXTRAIRE'  evu 'ORDO' ;
  728.  evu  = 'EVOL' 'MANU' 'x' lxx 'u' lu ;
  729.  tu   = 'CHAINE' 'u ' titolo ;
  730.  
  731. * v
  732.  
  733.  evv  = 'EVOL' 'CHPO' VNT LIGEVO ;
  734.  lv   = 'EXTRAIRE'  evv 'ORDO' ;
  735.  evv  = 'EVOL' 'MANU' 'x' lxx 'v' lv ;
  736.  tv   = 'CHAINE' 'v ' titolo ;
  737.  
  738.  
  739. * p
  740.  
  741.  evp  = 'EVOL' 'CHPO' PN LIGEVO ;
  742.  lp   = 'EXTRAIRE'  evp 'ORDO' ;
  743.  evp  = 'EVOL' 'MANU' 'x' lxx 'p' lp ;
  744.  tp   = 'CHAINE' 'p  ' titolo ;
  745.  
  746.  
  747. * T
  748.  
  749.  evT  = 'EVOL' 'CHPO' TN LIGEVO ;
  750.  lT   = 'EXTRAIRE'  evT 'ORDO' ;
  751.  evT  = 'EVOL' 'MANU' 'x' lxx 'T' lT ;
  752.  tT   = 'CHAINE' 'T  ' titolo ;
  753.  
  754.  'SI' GRAPH ;
  755.     'DESSIN' evro 'TITRE' tro 'XBOR' x0 x1 'MIMA' ;
  756.     'DESSIN' evp  'TITRE' tp 'XBOR' x0 x1  'MIMA' ;
  757.     'DESSIN' evu 'TITRE'  tu 'XBOR' x0 x1  'MIMA' ;
  758.     'DESSIN' evv 'TITRE'  tv 'XBOR' x0 x1  'MIMA' ;
  759.     'DESSIN' evT  'TITRE' tT 'XBOR' x0 x1  'MIMA' ;
  760.  'FINSI' ;
  761.  
  762.  
  763. * Test
  764.  
  765.  PVN = 3.2E6 ;
  766.  PCJ = 1.1E6 ;
  767.  RHOVN = 3.7 ;
  768.  RHOCJ = 1.4 ;
  769.  
  770.  
  771.  PMAX = 'MAXIMUM' PN ;
  772.  'SI' ((PMAX '>' PVN) 'OU' (PMAX '<' PCJ)) ;
  773.     'MESSAGE' 'Pression = ???' ;
  774.     'ERREUR' 5 ;
  775.  'FINSI' ;
  776.  
  777.  RHOMAX = 'MAXIMUM' RN ;
  778.  'SI' ((RHOMAX '>' RHOVN) 'OU' (RHOMAX '<' (0.85 '*' RHOCJ))) ;
  779.     'MESSAGE' 'Densite = ???' ;
  780.     'ERREUR' 5 ;
  781.  'FINSI' ;
  782.  
  783.  'FIN' ;
  784.  
  785.  
  786.  
  787.  
  788.  
  789.  
  790.  
  791.  
  792.  
  793.  
  794.  
  795.  

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