Télécharger la plaquette Cast3M


fr en

Télécharger konv_resi_dem3D_stationaryshock_21.dgibi

Retour à la liste

Numérotation des lignes : 

  1. * fichier :  konv_resi_dem3D_stationaryshock_21.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4. ***********************************************************
  5. **** APPROCHE VF "Cell-Centred Formulation" pour la    ****
  6. **** solution des                                      **** 
  7. **** Equations d'Euler pour un gaz parfait             ****
  8. **** Approche DEM pour la combustion                   ****
  9. **** OPERATEURS PRET, KONV                             ****
  10. ****                                                   ****
  11. **** 3D                                                ****
  12. ****                                                   ****
  13. **** Stationary shock  (phase 1 '-' vs phase 2)        ****
  14. ****                                                   ****
  15. **** Methodes: VLH, SS                                 ****
  16. ****                                                   ****
  17. **** A. BECCANTINI DM2S/SFME/LTMF    DECEMBRE 2010     ****
  18. ***********************************************************
  19.  
  20.  'OPTION'  'DIME' 3 'ELEM' 'CUB8' 'ECHO' 0 'TRAC' 'X' ;
  21.  
  22. *
  23. *** GRAPH
  24. *
  25.  
  26.  GRAPH = FAUX ;
  27. * GRAPH = VRAI ;
  28.  
  29. **** Cas homogene : la table proprieté de gaz
  30.  
  31.  
  32.  
  33. *************************************************
  34. **** The table for the properties of the gas ****
  35. *************************************************
  36. *
  37.  PGAS = 'TABLE' ;
  38. *
  39. **** Order of the polynomial order for cv = cv(T) 
  40. *    For T > TMAX, cv(T) = cv(Tmax)
  41. * 
  42.  PGAS . 'TMAX' = 6000.0 ;
  43.  PGAS . 'NORD' = 4 ;
  44. *
  45. **** Species involved in the mixture (before or after
  46. *    the chemical reaction)
  47. *
  48.  PGAS . 'SPECIES' = 'MOTS' 'H2  ' 'O2  ' 'H2O ' 'N2  ' ;
  49. *
  50. *
  51. **** Coefficient of the chemical reaction.
  52. *    Note that for the first species this coefficient should be positive
  53. *    Normal, we take it equal to 1.
  54. *
  55. *    H2 '+' 0.5 O2 ---> H2O
  56. *
  57.  PGAS . 'CHEMCOEF' = 'PROG' 1.0 0.5 -1.0 0.0 ;
  58. * 
  59. **** Mass fraction of the first species before and after the combustion
  60. *    Final mass fractions of the species with positive coefficients.
  61. *    Final mass fractions of the species with non-positive coefficient.
  62. *    The mass fraction of the last species is not given.
  63. *
  64.  PGAS . 'MASSFRA' = 'PROG' 0.285219E-01  0.964039E-11  0.765104E-10
  65.                            0.127442E-10 ;
  66. *
  67. **** Coef with the gas properties
  68. *
  69.  PGAS .  'H2  ' = 'TABLE'  ;
  70.  PGAS .  'H2O ' = 'TABLE'  ;
  71.  PGAS .  'N2  ' = 'TABLE'  ;
  72.  PGAS .  'O2  ' = 'TABLE'  ;
  73. *
  74. **** Runiv (J/mole/K)
  75. *
  76.  PGAS . 'RUNIV' = 8.31441 ;
  77. *
  78. **** W (kg/mole). Gas constant (J/kg/K = Runiv/W)
  79. *
  80.  PGAS .  'H2  ' . 'W' =  2. * 1.00797E-3 ;
  81.  PGAS .  'O2  ' . 'W' =  2. * 15.9994E-3 ;
  82.  PGAS .  'H2O ' . 'W' = (PGAS .  'H2  ' . 'W' ) '+'
  83.     (0.5 * (PGAS .  'O2  ' . 'W' )) ;
  84.  PGAS .  'N2  ' . 'W' = 2 * 14.0067E-3 ;
  85. *
  86. **** Polynomial coefficients 
  87. *
  88.  PGAS .  'H2  ' . 'A' = 'PROG'  9834.91866 0.54273926 0.000862203836
  89.                                -2.37281455E-07 1.84701105E-11 ;
  90.  PGAS .  'H2O ' . 'A' = 'PROG' 1155.95625 0.768331151 -5.73129958E-05
  91.                               -1.82753232E-08 2.44485692E-12 ;
  92.  PGAS .  'N2  ' . 'A' = 'PROG' 652.940766 0.288239099 -7.80442298E-05
  93.                               8.78233606E-09 -3.05514485E-13 ;
  94.  PGAS .  'O2  ' . 'A' = 'PROG' 575.012333  0.350522002 -0.000128294865
  95.                               2.33636971E-08 -1.53304905E-12;
  96. *
  97. **** Formation enthalpies (energies) at 0K (J/Kg)
  98. *
  99.  PGAS .  'H2  ' . 'H0K' = -4.195D6 ;
  100.  PGAS .  'H2O ' . 'H0K' = -1.395D7 ;
  101.  PGAS .  'N2  ' . 'H0K' = -2.953D5 ;
  102.  PGAS .  'O2  ' . 'H0K' = -2.634D5 ;
  103.  
  104.  
  105. *************************** 
  106. ***** DOMAINE SPATIAL  ****
  107. ***************************
  108.  
  109.  A1 = 0.0D0 0.0D0 0.0;
  110.  A2 = 1.0D0 0.0D0 0.0;
  111.  A3 = 2.0D0 0.0D0 0.0;
  112.  A4 = 2.0D0 1.0D0 0.0;
  113.  A5 = 1.0D0 1.0D0 0.0;
  114.  A6 = 0.0D0 1.0D0 0.0;
  115.  
  116.  L12 = A1 'DROIT' 1 A2;
  117.  L23 = A2 'DROIT' 1 A3;
  118.  L34 = A3 'DROIT' 1 A4;
  119.  L45 = A4 'DROIT' 1 A5;
  120.  L56 = A5 'DROIT' 1 A6;
  121.  L61 = A6 'DROIT' 1 A1;
  122.  L25 = A2 'DROIT' 1 A5;
  123.  
  124.  
  125.  BAS1 = 'DALL' L12 L25 L56 L61
  126.          'PLANE';
  127.  DOM1 = 'VOLUME' BAS1 'TRANSLATION' 1 (0.0 0.0 1.0) ;
  128.  BAS2 = 'DALL' L23 L34 L45 ('INVERSE' L25)
  129.          'PLANE';
  130.  DOM2 = 'VOLUME' BAS2 'TRANSLATION' 1 (0.0 0.0 1.0) ;
  131.  
  132.  DOM10 = DOM1 ;
  133.  DOM20 = DOM2 ;
  134.  DOMTOT = DOM1 'ET' DOM2 ;
  135.  'ELIMINATION' 0.0001 DOMTOT ;
  136.  
  137. *
  138. *** Point ou on controlle la consistence
  139. *
  140.  
  141.  P10 = 1.0 0.5 0.5 ;
  142.  
  143.  $DOMTOT = 'MODELISER' DOMTOT 'EULER';
  144.  $DOM1 = 'MODELISER'  DOM1 'EULER';
  145.  $DOM2 = 'MODELISER'  DOM2 'EULER';
  146.  TDOMTOT = 'DOMA'  $DOMTOT 'VF';
  147.  TDOM1 = 'DOMA'  $DOM1 'VF';
  148.  TDOM2 = 'DOMA'  $DOM2 'VF';
  149.  MDOM1 = TDOM1 . 'QUAF' ;
  150.  MDOM2 = TDOM2 . 'QUAF' ;
  151.  MDOMTOT = TDOMTOT . 'QUAF' ;
  152.  
  153.  'ELIM' (MDOMTOT 'ET' MDOM1 'ET' MDOM2) 1.E-6 ;
  154.  
  155.  
  156.  'SI' GRAPH;
  157.     'TRACER' (('DOMA' $DOMTOT 'MAILLAGE') 'ET'  
  158.        ('DOMA' $DOMTOT 'FACEL') 'ET' P1) 'TITRE' 'Domaine et FACEL';
  159.  'FINSI' ;
  160.  
  161. *
  162.  EPS = 1.0D-16 ;
  163.  K0 = 45.2 ;
  164. * 
  165. *******************************************
  166. *** Etats gauche et droite ****************
  167. *******************************************
  168. *
  169. * There are 3 fluxes to control :
  170. * 1 -> 2  can be checked by considering an exact reactive shock
  171. * 2 -> 1  can be checked by considering an exact reactive shock
  172. * 1 -> 1  can be checked via consistency (or via an exact non-reactive
  173. *         or by comparing the results with a non-reactive one)
  174. * 2 -> 2  can be checked via consistency (or via an exact non-reactive
  175. *         or by comparing the results with a non-reactive one)
  176. *
  177. * Case 1 -> 1
  178. * Constant state
  179. *
  180.  epsv = 0.0 ;
  181.  tg1 = 327.0 ;
  182.  td1 = 327.0 ;
  183.  pg1 = 139000.0 ;
  184.  pd1 = 139000.0 ;
  185.  ung1 = -45.2 '-' epsv ;
  186.  und1 = -45.2 '-' epsv ;
  187.  utg1 = 0.0 ;
  188.  utd1 = 0.0 ;
  189.  uvg1 = 0.0 ;
  190.  uvd1 = 0.0 ;
  191. * utg1 = 0.0 ;
  192. * utd1 = 0.0 ;
  193.  tg2 = 2518.582685186033 ;
  194.  td2 = 2518.582685186033 ;
  195.  pg2 = 125279.9833247705 ;
  196.  pd2 = 125279.9833247705 ;
  197.  ung2 = -329.121088555842 '-' epsv ;
  198.  und2 = -329.121088555842 '-' epsv ;
  199.  utg2 = 0.0 ;
  200.  utd2 = 0.0 ;
  201.  uvg2 = 0.0 ;
  202.  uvd2 = 0.0 ;
  203.  
  204.  TN1 = 'PROG' tg1 td1 ;
  205.  TN2 = 'PROG' tg2 td2 ; 
  206.  PN1 = 'PROG' pg1 pd1 ;
  207.  PN2 = 'PROG' pg2 pd2 ;
  208.  UN1 = 'PROG' ung1 und1 ;
  209.  UN2 = 'PROG' ung2 und2 ;
  210.  UT1 = 'PROG' utg1 utd1 ;
  211.  UT2 = 'PROG' utg2 utd2 ;
  212.  UV1 = 'PROG' uvg1 uvd1 ;
  213.  UV2 = 'PROG' uvg2 uvd2 ;
  214.  
  215.  ALPHA1 = 'PROG' 0.0 1.0 ;
  216.  ALPHA2 = ('PROG' ('DIME' ALPHA1) '*' 1.0) '-' ALPHA1 ;
  217.  
  218. * 
  219. *
  220. * ETHER = int_0^T cv(T') dT'           T <  TMAX
  221. *       = int_0^TMAX cv(T') dT' '+'
  222. *         cv(TMAX)                     T >= TMAX
  223. *
  224.  
  225.  ESP1 = 'EXTRAIRE' (PGAS . 'SPECIES') 1 ;
  226. * DY1 = y_i - y_f for the species 1 
  227.  DY1 = (('EXTRAIRE' (PGAS . 'MASSFRA') 1) '-'
  228.    ('EXTRAIRE' (PGAS . 'MASSFRA') 2)) ;
  229.  COEF1 = ('EXTRAIRE' (PGAS . 'CHEMCOEF') 1) '*'
  230.     (PGAS . ESP1 . 'W') ;
  231.  YFINPH1 = 1.0 ;
  232.  YFINPH2 = 1.0 ;
  233.  'SI' (COEF1 > 0) ;
  234.     YPH2 = 'EXTRAIRE' (PGAS . 'MASSFRA') 2 ;
  235.     YPH1 = YPH2 '+' DY1 ;
  236.  'SINON' ;
  237.     YPH1 = 'EXTRAIRE' (PGAS . 'MASSFRA') 2 ;
  238.     YPH2 = YPH1 '-' DY1 ;
  239.  'FINSI' ;
  240.  YFINPH1 = YFINPH1 '-' YPH1 ;
  241.  YFINPH2 = YFINPH2 '-' YPH2 ;
  242.  PRYPH1 = 'PROG' YPH1 ;
  243.  PRYPH2 = 'PROG' YPH2 ;
  244.  'REPETER' BLESP (('DIME' (PGAS . 'SPECIES')) '-' 2) ;
  245.     ESP = 'EXTRAIRE' (PGAS . 'SPECIES') (&BLESP '+' 1)  ;
  246.     COEF = ('EXTRAIRE' (PGAS . 'CHEMCOEF') (&BLESP '+' 1)) 
  247.        '*' (PGAS . ESP . 'W') ;
  248.     DY = (DY1 * (COEF '/' COEF1)) ;
  249.     'SI' (COEF > 0) ;
  250.        YPH2 = 'EXTRAIRE' (PGAS . 'MASSFRA') (&BLESP '+' 2) ;
  251.        YPH1 = YPH2 '+' DY ;
  252.     'SINON' ;
  253.        YPH1 = 'EXTRAIRE' (PGAS . 'MASSFRA') (&BLESP '+' 2) ;
  254.        YPH2 = YPH1 '-' DY ;
  255.     'FINSI' ;
  256.     PRYPH1 = PRYPH1 'ET' ('PROG' YPH1) ;
  257.     PRYPH2 = PRYPH2 'ET' ('PROG' YPH2) ;
  258.     YFINPH1 = YFINPH1 '-' YPH1 ;
  259.     YFINPH2 = YFINPH2 '-' YPH2 ;
  260.  'FIN' BLESP ;
  261.  PRYPH1 = PRYPH1 'ET' ('PROG' YFINPH1) ; 
  262.  PRYPH2 = PRYPH2 'ET' ('PROG' YFINPH2) ;
  263.  'LISTE' PRYPH1 ;
  264.  'LISTE' PRYPH2 ; 
  265. *
  266.  TMAX = ('PROG' ('DIME' ALPHA1) '*' (PGAS . 'TMAX')) ;
  267. * TCAL1 = MIN TN1, TMAX 
  268.  TCAL1 = 0.5D0 '*' ((TMAX '+' TN1) '-' ('ABS' (TN1 '-' TMAX))) ;
  269.  DTN1 = TN1 '-' TCAL1 ;
  270. * TCAL1 = MIN TN1, TMAX 
  271.  TCAL2 = 0.5D0 '*' ((TMAX '+' TN2) '-' ('ABS' (TN2 '-' TMAX))) ;
  272.  DTN2 = TN2 '-' TCAL2 ;
  273. *
  274. * Internal energy (J/kg in SI)
  275. *
  276.  ETHER1 = ('PROG' ('DIME' ALPHA1) '*' 0.0) ;
  277.  CV1 = ('PROG' ('DIME' ALPHA1) '*' 0.0) ;
  278.  ETHER2 = ('PROG' ('DIME' ALPHA1) '*' 0.0) ;
  279.  CV2 = ('PROG' ('DIME' ALPHA1) '*' 0.0) ;
  280.  FUNTN1 = ('PROG' ('DIME' ALPHA1) '*' 1.0) ;
  281.  FUNTN2 = ('PROG' ('DIME' ALPHA1) '*' 1.0) ;
  282.  'REPETER' BLPO ((PGAS . 'NORD') '+' 1) ;
  283.     'REPETER' BLESP ('DIME' (PGAS . 'SPECIES')) ;
  284.        ESP = 'EXTRAIRE' (PGAS . 'SPECIES') &BLESP ;
  285.        YCEL1 = 'EXTRAIRE' PRYPH1 &BLESP ;
  286.        YCEL2 = 'EXTRAIRE' PRYPH2 &BLESP ;
  287.        AA = 'EXTRAIRE' (PGAS . ESP . 'A') &BLPO ;
  288.        DCV1 = (AA * YCEL1 * FUNTN1) ;
  289.        DCV2 = (AA * YCEL2 * FUNTN2) ;
  290.        CV1 = CV1 '+' DCV1 ;
  291.        CV2 = CV2 '+' DCV2 ;
  292.        ETHER1 = ETHER1 '+' (DCV1 * TCAL1 '/' (&BLPO)) ;
  293.        ETHER2 = ETHER2 '+' (DCV2 * TCAL2 '/' (&BLPO)) ;
  294.     'FIN' BLESP ;
  295.     FUNTN1 = FUNTN1 '*' TCAL1 ;
  296.     FUNTN2 = FUNTN2 '*' TCAL2 ;
  297.  'FIN' BLPO ;
  298.  ETHER1 = ETHER1 '+' (CV1 '*' DTN1) ;
  299.  ETHER2 = ETHER2 '+' (CV2 '*' DTN2) ;
  300. *
  301. * Formation energy/enthalpy (J/kg in SI) and gas constant (J/kg/K)
  302. *
  303.  EFORM1 = 0.0 ;
  304.  EFORM2 = 0.0 ;
  305.  RGAS1 = 0.0 ;
  306.  RGAS2 = 0.0 ;
  307.  'REPETER' BLESP ('DIME' (PGAS . 'SPECIES')) ;
  308.     ESP = 'EXTRAIRE' (PGAS . 'SPECIES') &BLESP ;
  309.     YCEL1 = 'EXTRAIRE' PRYPH1 &BLESP ;
  310.     YCEL2 = 'EXTRAIRE' PRYPH2 &BLESP ;
  311.     EFORM1 = EFORM1 '+' (YCEL1 * (PGAS . ESP . 'H0K')) ;
  312.     EFORM2 = EFORM2 '+' (YCEL2 * (PGAS . ESP . 'H0K')) ;
  313.     RGAS1 = RGAS1 '+' (YCEL1 *  (PGAS . 'RUNIV') '/' 
  314.       (PGAS . ESP . 'W')) ;
  315.     RGAS2 = RGAS2 '+' (YCEL2 *  (PGAS . 'RUNIV') '/' 
  316.       (PGAS . ESP . 'W')) ;
  317.  'FIN' BLESP ;
  318.  EFORM1 = 'PROG' ('DIME' ETHER1) '*' EFORM1 ;
  319.  EFORM2 = 'PROG' ('DIME' ETHER2) '*' EFORM2 ;
  320.  'LISTE' rgas1 ;
  321.  'LISTE' rgas2 ;
  322.  
  323. *
  324. * Computation of the conservative variables
  325. *
  326.  RN1 = PN1 '/' (RGAS1 '*' TN1) ;
  327.  RN2 = PN2 '/' (RGAS2 '*' TN2) ;
  328.  GNX1 = RN1 * UN1 ;
  329.  GNY1 = RN1 * UT1 ;
  330.  GNX2 = RN2 * UN2 ;
  331.  GNY2 = RN2 * UT2 ;
  332.  ECIN1 = 0.5D0 '*' ('PROG'
  333.     ((ung1 * ung1) '+' (utg1 * utg1) '+' (uvg1 * uvg1))
  334.     ((und1 * und1) '+' (utd1 * utd1) '+' (uvd1 * uvd1))
  335.     );
  336.  ECIN2 = 0.5D0 '*' ('PROG'
  337.     ((ung2 * ung2) '+' (utg2 * utg2) '+' (uvg2 * uvg2))
  338.     ((und2 * und2) '+' (utd2 * utd2) '+' (uvd2 * uvd2)));
  339.  RETN1 = RN1 '*' (ETHER1 '+' ECIN1 '+' EFORM1) ; 
  340.  RETN2 = RN2 '*' (ETHER2 '+' ECIN2 '+' EFORM2) ;
  341. *
  342. * Computation of the resi contribution in (n,t)
  343. *
  344. * State 1
  345.  retg1 = MAXI RETN1 ;
  346.  rhog1 = 'MAXIMUM' RN1 ;
  347.  fA11gd = ung1 '*' rhog1 ;
  348.  fA12gd = (fA11gd '*' ung1) '+' pg1 ;
  349.  fA13gd = fA11gd '*' utg1 ;
  350.  fA14gd = ung1 '*' (retg1 '+' pg1);
  351. *
  352. * State 2 
  353.  retg2 = MAXI RETN2 ;
  354.  rhog2 = 'MAXIMUM' RN2 ;
  355.  fA21gd = ung2 '*' rhog2 ;
  356.  fA22gd = (fA21gd '*' ung2) '+' pg2 ;
  357.  fA23gd = fA21gd '*' utg2 ;
  358.  fA24gd = ung2 '*' (retg2 '+' pg2);
  359. *
  360.  
  361. ****************************************************
  362. ****************************************************
  363. ******** Boucle sur les angles             *********
  364. ****************************************************
  365. ****************************************************
  366.  
  367.  DANGLE = 360 '/' 7.15;
  368.  ANGLE = 0.0  ;
  369.  
  370.  'REPETER' BLOC 8;
  371.  
  372. *
  373. *** Rotation
  374. *
  375.  
  376.     ANGLE = ANGLE '+' DANGLE;
  377.     ORIG = 0.0D0 0.0D0 0.0 ;
  378.     P1   = 1.0 1.0 1.0 ;
  379.  
  380.     'MESSAGE' ;
  381.     'MESSAGE'  (CHAIN 'Angle de rotation= ' ANGLE);
  382.     'MESSAGE' ;
  383.  
  384.     DOM1   = DOM10   'TOURNER' ANGLE ORIG P1 ;
  385.     DOM2   = DOM20   'TOURNER' ANGLE ORIG P1 ;
  386.     P1FAC  = P10     'TOURNER' ANGLE ORIG P1 ;
  387.  
  388.     DOMTOT = DOM1 ET DOM2;
  389.     'ELIMINATION' DOMTOT 1D-6;
  390.  
  391.     $DOMTOT = 'MODELISER' DOMTOT 'EULER';
  392.     $DOM1 = 'MODELISER'  DOM1 'EULER';
  393.     $DOM2 = 'MODELISER'  DOM2 'EULER';
  394.     TDOMTOT = 'DOMA'  $DOMTOT 'VF';
  395.     TDOM1 = 'DOMA'  $DOM1 'VF';
  396.     TDOM2 = 'DOMA'  $DOM2 'VF';
  397.     MDOM1 = TDOM1 . 'QUAF' ;
  398.     MDOM2 = TDOM2 . 'QUAF' ;
  399.     MDOMTOT = TDOMTOT . 'QUAF' ;
  400.  
  401.     'ELIM' (MDOMTOT 'ET' MDOM1 'ET' MDOM2) 1.E-6 ;
  402.     GRALP1 = 'NOMC' ('MOTS' 'UX' 'UY' 'UZ') (TDOMTOT . 'XXNORMAF')
  403.        ('MOTS' 'P1DX' 'P1DY' 'P1DZ') ;
  404.  
  405. *
  406.     VINF = 'MANUEL' 'CHPO' (TDOMTOT . 'CENTRE') 1 'SCAL' 100. ;
  407. *
  408. **** Redefinition de P1FAC dans $DOMTOT 'FACE'
  409. *
  410.  
  411.     P1FAC = ('DOMA' $DOMTOT 'FACE') 'POIN' 'PROC' P1FAC ;
  412.     XSURF = 'EXTRAIRE' ('DOMA' $DOMTOT 'XXSURFAC') P1FAC 'SCAL' ;
  413.     NX    = 'EXTRAIRE' (TDOMTOT . 'MATROT') P1FAC 'MX' ;
  414.     NY    = 'EXTRAIRE' (TDOMTOT . 'MATROT') P1FAC 'MY' ;
  415.     NZ    = 'EXTRAIRE' (TDOMTOT . 'MATROT') P1FAC 'MZ' ;
  416.     TX    = 'EXTRAIRE' (TDOMTOT . 'MATROT') P1FAC 'UX' ;
  417.     TY    = 'EXTRAIRE' (TDOMTOT . 'MATROT') P1FAC 'UY' ; 
  418.     TZ    = 'EXTRAIRE' (TDOMTOT . 'MATROT') P1FAC 'UZ' ;
  419.     VX    = 'EXTRAIRE' (TDOMTOT . 'MATROT') P1FAC 'RX' ;
  420.     VY    = 'EXTRAIRE' (TDOMTOT . 'MATROT') P1FAC 'RY' ;
  421.     VZ    = 'EXTRAIRE' (TDOMTOT . 'MATROT') P1FAC 'RZ' ;
  422.  
  423.     GEOPC    = ('DOMA' $DOMTOT 'FACEL') 'ELEM' 'APPUYE' 'LARGEMENT'
  424.         P1FAC ;
  425.     AA = 'CHANGER' GEOPC 'POI1' ;
  426.     PC1 = 'POIN' 1 AA ;
  427.     P2 = 'POIN' 1 ('DOMA' $DOMTOT 'CENTRE') ;
  428.     PC2 = 'POIN' 3 AA ;
  429. *   Si P3 n'existe pas, probleme en FACEL.     
  430.     'SI' ('NEG' PC1 P2) ;
  431.        ORIENT = -1 ;
  432.        PCD = PC1 ;
  433.        PCG = PC2 ;
  434.        XVOLD = 'MAXIMUM' (TDOM1 . 'XXVOLUM') ;
  435.        XVOLG = 'MAXIMUM' (TDOM2 . 'XXVOLUM') ;
  436.     'SINON' ;
  437.        ORIENT = 1 ;
  438.        PCD = PC2 ;
  439.        PCG = PC1 ;
  440.        XVOLG = 'MAXIMUM' (TDOM1 . 'XXVOLUM') ;
  441.        XVOLD = 'MAXIMUM' (TDOM2 . 'XXVOLUM') ;
  442.     'FINSI' ;
  443.  
  444.     uxg1 = ((ung1 '*' NX) '+' (utg1 '*' TX) '+' (uvg1 '*' VX))
  445.        '*' ORIENT ;
  446.     uyg1 = ((ung1 '*' NY) '+' (utg1 '*' TY) '+' (uvg1 '*' VY))
  447.        '*' ORIENT ;
  448.     uzg1 = ((ung1 '*' NZ) '+' (utg1 '*' TZ) '+' (uvg1 '*' VZ))
  449.        '*' ORIENT ;
  450.     uxg2 = ((ung2 '*' NX) '+' (utg2 '*' TX) '+' (uvg2 '*' VX))
  451.        '*' ORIENT ;
  452.     uyg2 = ((ung2 '*' NY) '+' (utg2 '*' TY) '+' (uvg2 '*' VY))
  453.        '*' ORIENT ;
  454.     uzg2 = ((ung2 '*' NZ) '+' (utg2 '*' TZ) '+' (uvg2 '*' VZ))
  455.        '*' ORIENT ;
  456.     uxd1 = ((und1 '*' NX) '+' (utd1 '*' TX) '+' (uvd1 '*' VX))
  457.        '*' ORIENT ;
  458.     uyd1 = ((und1 '*' NY) '+' (utd1 '*' TY) '+' (uvd1 '*' VY))
  459.        '*' ORIENT ;
  460.     uzd1 = ((und1 '*' NZ) '+' (utd1 '*' TZ) '+' (uvd1 '*' VZ))
  461.        '*' ORIENT ;
  462.     uxd2 = ((und2 '*' NX) '+' (utd2 '*' TX) '+' (uvd2 '*' VX))
  463.        '*' ORIENT ;
  464.     uyd2 = ((und2 '*' NY) '+' (utd2 '*' TY) '+' (uvd2 '*' VY))
  465.        '*' ORIENT ;
  466.     uzd2 = ((und2 '*' NZ) '+' (utd2 '*' TZ) '+' (uvd2 '*' VZ))
  467.        '*' ORIENT ;
  468. *    'OPTION' DONN 5 ;
  469.  
  470.     'SI' GRAPH;
  471.        'TRACER' (('DOMA' $DOMTOT 'MAILLAGE') 'ET'  
  472.           ('COULEUR' ('DOMA' $DOMTOT 'FACEL') 'ROUGE') 'ET' P1FAC)
  473.           'TITRE' 'Domaine et FACEL';
  474.     'FINSI' ;
  475.  
  476.     MAILIM = 'DIFF' ('DOMA' $DOMTOT 'FACE') ('MANUEL' 'POI1' P1FAC) ;
  477.  
  478. ***********************
  479. **** Les CHPOINTs  ****
  480. ***********************
  481.  
  482.     CHRN1   = ('MANUEL' 'CHPO' (TDOMTOT . 'CENTRE') 1 'SCAL' RN1) ;
  483.     CHRN2   = ('MANUEL' 'CHPO' (TDOMTOT . 'CENTRE') 1 'SCAL' RN2) ;
  484.  
  485.     CHVN1   = ('MANUEL' 'CHPO' (TDOM1 . 'CENTRE') 3 'UX' uxg1
  486.        'UY' uyg1 'UZ' uzg1) '+' 
  487.               ('MANUEL' 'CHPO' (TDOM2 . 'CENTRE') 3 'UX' uxd1
  488.        'UY' uyd1 'UZ' uzd1) ;
  489.     CHVN2   = ('MANUEL' 'CHPO' (TDOM1 . 'CENTRE') 3 'UX' uxg2
  490.        'UY' uyg2 'UZ' uzg2) '+' 
  491.               ('MANUEL' 'CHPO' (TDOM2 . 'CENTRE') 3 'UX' uxd2
  492.        'UY' uyd2 'UZ' uzd2) ;
  493.  
  494.     CHGN1 = CHRN1 '*' CHVN1 ;
  495.     CHGN2 = CHRN2 '*' CHVN2 ;
  496.  
  497.     CHRET1   = ('MANUEL' 'CHPO' (TDOMTOT . 'CENTRE') 1 'SCAL' RETN1) ;
  498.     CHRET2   = ('MANUEL' 'CHPO' (TDOMTOT . 'CENTRE') 1 'SCAL' RETN2) ;
  499.  
  500.     CHTN1   = ('MANUEL' 'CHPO' (TDOMTOT . 'CENTRE') 1 'SCAL' TN1) ;
  501.     CHTN2   = ('MANUEL' 'CHPO' (TDOMTOT . 'CENTRE') 1 'SCAL' TN2) ;
  502.  
  503.     CHAL1 = ('MANUEL' 'CHPO' (TDOMTOT . 'CENTRE') 1 'SCAL' ALPHA1) ;
  504.     CHAL2 = ('MANUEL' 'CHPO' (TDOMTOT . 'CENTRE') 1 'SCAL' ALPHA2) ;
  505.  
  506.     R1 R2 V1 V2 P1 P2 T1 T2 = 'PRIM' 'DEM' PGAS CHAL1 CHAL2
  507.       (CHAL1 * CHRN1) (CHAL2 * CHRN2) (CHAL1 * CHGN1)
  508.       (CHAL2 * CHGN2) (CHAL1 * CHRET1) (CHAL2 * CHRET2)
  509.        CHTN1 CHTN2 EPS ;
  510.  
  511.     GRADAL1 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 3
  512.       'P1DX' 0.0 'P1DY' 0.0 'P1DZ' 0.0 ;
  513.     GRADAL2 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 3
  514.       'P1DX' 0.0 'P1DY' 0.0 'P1DZ' 0.0 ;
  515.     LIMAL1 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 1
  516.       'P1' 0.0 ; 
  517.     LIMAL2 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 1
  518.       'P1' 0.0 ; 
  519.     GRADR1 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 3
  520.       'P1DX' 0.0 'P1DY' 0.0 'P1DZ' 0.0;
  521.     GRADR2 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 3
  522.       'P1DX' 0.0 'P1DY' 0.0 'P1DZ' 0.0;
  523.     LIMR1 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 1
  524.       'P1' 0.0 ; 
  525.     LIMR2 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 1
  526.       'P1' 0.0 ; 
  527.     GRADV1 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 9
  528.       'P1DX' 0.0 'P1DY' 0.0 'P1DZ' 0.0
  529.       'P2DX' 0.0 'P2DY' 0.0 'P2DZ' 0.0
  530.       'P3DX' 0.0 'P3DY' 0.0 'P3DZ' 0.0 ;
  531.     GRADV2 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 9
  532.       'P1DX' 0.0 'P1DY' 0.0 'P1DZ' 0.0
  533.       'P2DX' 0.0 'P2DY' 0.0 'P2DZ' 0.0
  534.       'P3DX' 0.0 'P3DY' 0.0 'P3DZ' 0.0 ;
  535.     LIMV1 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 3
  536.       'P1' 0.0 'P2' 0.0 'P3' 0.0 ; 
  537.     LIMV2 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 3
  538.       'P1' 0.0 'P2' 0.0 'P3' 0.0 ; 
  539.     GRADP1 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 3
  540.       'P1DX' 0.0 'P1DY' 0.0 'P1DZ' 0.0;
  541.     GRADP2 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 3
  542.       'P1DX' 0.0 'P1DY' 0.0 'P1DZ' 0.0;
  543.     LIMP1 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 1
  544.       'P1' 0.0 ; 
  545.     LIMP2 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 1
  546.       'P1' 0.0 ; 
  547.  
  548. *
  549. ****  L'operateur  'PRET'
  550. *
  551.  
  552.     CHFAL1 CHFAL2 CHFRN1 CHFRN2 CHFVN1 CHFVN2 CHFPN1 CHFPN2 =
  553.      'PRET' 'DEM'  $DOMTOT
  554.       CHAL1  (0.0 * GRADAL1)  LIMAL1
  555.       CHAL2  (0.0 * GRADAL2)  LIMAL2
  556.       R1     (0.0 * GRADR1)   LIMR1
  557.       R2     (0.0 * GRADR2)   LIMR2
  558.       V1     (0.0 * GRADV1)   LIMV1
  559.       V2     (0.0 * GRADV2)   LIMV2
  560.       P1     (0.0 * GRADP1)   LIMP1
  561.       P2     (0.0 * GRADP2)   LIMP2 ;
  562. *
  563.     SI FAUX ;
  564.        CHFRO1  = 'REDU' CHFRN1 geopc ;
  565.        CHFV1   = 'REDU' CHFVN1 geopc ;
  566.        CHFP1   = 'REDU' CHFPN1 geopc ;
  567.        CHFRO2  = 'REDU' CHFRN2 geopc ;
  568.        CHFV2   = 'REDU' CHFVN2 geopc ;
  569.        CHFP2   = 'REDU' CHFPN2 geopc ;
  570.        'LISTE' CHFRO1 ;
  571.        'LISTE' CHFV1 ;
  572.        'LISTE' CHFP1 ;
  573.        'LISTE' CHFRO2 ;
  574.        'LISTE' CHFV2 ;
  575.        'LISTE' CHFP2 ;
  576.        'OPTION' donn 5 ;
  577.     'FINSI' ;
  578.  
  579.  
  580. *
  581. **** L'operateur 'KONV'
  582. *
  583.  
  584.     'REPETER' BLMETO 3 ;
  585.  
  586.        'SI' ('EGA' &BLMETO 1) ;
  587.            METO = 'SS'  ;
  588.        'FINSI' ;
  589.        'SI' ('EGA' &BLMETO 2) ;
  590.            METO = 'VLH' ;
  591.        'FINSI' ;
  592.        'SI' ('EGA' &BLMETO 3) ;
  593.            METO = 'AUSMPUP'  ;
  594.        'FINSI' ;
  595.  
  596.        'MESSAGE' ('CHAINE' 'METO = ' METO) ;
  597.        LISTINC1 = ('MOTS' 'ALF1' 'RN1' 'RNX1' 'RUY1' 'RUZ1' 'RET1'
  598.           'ALF2' 'RN2' 'RUX2' 'RUY2' 'RUZ2' 'RET2') ;
  599.  
  600.        SI ('EGA' METO 'AUSMPUP') ;
  601.           CHPRES DT SURF = 'KONV' 'VF' 'DEM' 'RESI' METO 'CONS'
  602.           $DOMTOT PGAS LISTINC1 CHAL1 CHAL2 CHFAL1 CHFAL2 CHFRN1 CHFRN2
  603.           CHFVN1 CHFVN2 CHFPN1 CHFPN2 K0 GRALP1 EPS MAILIM VINF VINF ;
  604.        'SINON' ;
  605.           CHPRES DT SURF = 'KONV' 'VF' 'DEM' 'RESI' METO 'CONS'
  606.           $DOMTOT PGAS LISTINC1 CHAL1 CHAL2 CHFAL1 CHFAL2 CHFRN1 CHFRN2
  607.           CHFVN1 CHFVN2 CHFPN1 CHFPN2 K0 GRALP1 EPS MAILIM ;
  608.        'FINSI' ;
  609.  
  610.  
  611. *       'LISTE' ('CHAINE' 'FORMAT' '(E16.10)' ('EXTRAIRE' ETHER1 1)) ;
  612. *       'LISTE' ('CHAINE' 'FORMAT' '(E16.10)' ('EXTRAIRE' ECIN1 1)) ;
  613. *       'LISTE' ('CHAINE' 'FORMAT' '(E16.10)' ('EXTRAIRE' EFORM1 1)) ;
  614. *       'LISTE' ('CHAINE' 'FORMAT' '(E16.10)' (pg1)) ;
  615. *       'LISTE' ('CHAINE' 'FORMAT' '(E16.10)' ('EXTRAIRE' RN1 1)) ;
  616. *       'LISTE' ('CHAINE' 'FORMAT' '(E16.10)' (ung1)) ;
  617. *       'OPTION' DONN 5 ;
  618.  
  619. *
  620. *      We verify that the residum contribution of P1FAC is the one
  621. *      of a reactive stationary shock.
  622. *      
  623.        RA11 = 'EXCO' ('EXTRAIRE' LISTINC1 1) CHPRES ;
  624.        RA12  = 'EXCO' ('EXTRAIRE' LISTINC1 2) CHPRES ;
  625.        RA13X = 'EXCO' ('EXTRAIRE' LISTINC1 3) CHPRES ;
  626.        RA13Y = 'EXCO' ('EXTRAIRE' LISTINC1 4) CHPRES ;
  627.        RA13Z = 'EXCO' ('EXTRAIRE' LISTINC1 5) CHPRES ;
  628.        RA13N = (RA13X '*' NX) '+' (RA13Y * NY) '+' (RA13Z * NZ) ;
  629.        RA13T = (RA13X '*' TX) '+' (RA13Y * TY) '+' (RA13Z * TZ) ;
  630.        RA13V = (RA13X '*' VX) '+' (RA13Y * VY) '+' (RA13Z * VZ) ;
  631.        RA14 = 'EXCO' ('EXTRAIRE' LISTINC1 6) CHPRES ;
  632. *
  633.        RA21 = 'EXCO' ('EXTRAIRE' LISTINC1 7) CHPRES ;
  634.        RA22  = 'EXCO' ('EXTRAIRE' LISTINC1 8) CHPRES ;
  635.        RA23X = 'EXCO' ('EXTRAIRE' LISTINC1 9) CHPRES ;
  636.        RA23Y = 'EXCO' ('EXTRAIRE' LISTINC1 10) CHPRES ;
  637.        RA23Z = 'EXCO' ('EXTRAIRE' LISTINC1 11) CHPRES ;
  638.        RA23N = (RA23X '*' NX) '+' (RA23Y * NY) '+' (RA23Z * NZ) ;
  639.        RA23T = (RA23X '*' TX) '+' (RA23Y * TY) '+' (RA23Z * TZ) ;
  640.        RA23V = (RA23X '*' VX) '+' (RA23Y * VY) '+' (RA23Z * VZ) ;
  641.        RA24 = 'EXCO' ('EXTRAIRE' LISTINC1 12) CHPRES ;
  642. *
  643.        ERRO = 1D-7  ;
  644.        HORROR = 'MAXIMUM' RA11 'ABS' ; 
  645.        LOGI1 = HORROR < ERRO ;
  646.        'SI' ('NON' LOGI1) ;
  647.           'MESSAGE' ('CHAINE' 'Erreur alpha ' HORROR ' > ' ERRO);
  648.        'FINSI' ;
  649.        HORROR = 'MAXIMUM' RA21 'ABS' ; 
  650.        LOGI2 = HORROR < ERRO ;
  651.        'SI' ('NON' LOGI2) ;
  652.           'MESSAGE' ('CHAINE' 'Erreur alpha ' HORROR ' > ' ERRO);
  653.        'FINSI' ;
  654.        LOGI1 = LOGI1 'ET' LOGI2;
  655. *
  656.        dRA12 = ('EXTRAIRE' RA12 'SCAL' PCD) '*' XVOLD ;
  657.        gRA12 = ('EXTRAIRE' RA12 'SCAL' PCG) '*' XVOLG ;
  658.        dRA22 = ('EXTRAIRE' RA22 'SCAL' PCD) '*' XVOLD ;
  659.        gRA22 = ('EXTRAIRE' RA22 'SCAL' PCG) '*' XVOLG ;
  660.        ERRO = 1.0D-8 * fA11gd 'ABS' ;
  661.        HORROR = ('MAXIMUM' ('PROG' gRA12 (dRA12 '-' fA11gd)
  662.           (gRA22 '+' fA21gd) dRA22) 'ABS') ;
  663.        LOGI2 = HORROR < ERRO ;
  664.        'SI' ('NON' LOGI2) ;
  665.           'MESSAGE' ('CHAINE' 'Erreur density ' HORROR ' > ' ERRO);
  666.        'FINSI' ;
  667.        LOGI1 = LOGI1 'ET' LOGI2;
  668. *
  669.        dRA13N = ('EXTRAIRE' RA13N 'SCAL' PCD) '*' XVOLD  ;
  670.        gRA13N = ('EXTRAIRE' RA13N 'SCAL' PCG) '*' XVOLG  ;
  671.        dRA23N = ('EXTRAIRE' RA23N 'SCAL' PCD) '*' XVOLD  ;
  672.        gRA23N = ('EXTRAIRE' RA23N 'SCAL' PCG) '*' XVOLG  ;
  673.        ERRO = 1.0D-8 * fA12gd 'ABS' ;
  674.        HORROR = ('MAXIMUM' ('PROG' gRA13N (dRA13N '-' fA12gd)
  675.           (gRA23N '+' fA22gd) dRA23N) 'ABS') ;
  676.        LOGI2 = HORROR < ERRO ;
  677.        'SI' ('NON' LOGI2) ;
  678.           'MESSAGE' ('CHAINE' 'Erreur qdm_n ' HORROR ' > ' ERRO);
  679.        'FINSI' ;
  680.        LOGI1 = LOGI1 'ET' LOGI2;
  681. *
  682.        dRA13T = ('EXTRAIRE' RA13T 'SCAL' PCD) '*' XVOLD  ;
  683.        gRA13T = ('EXTRAIRE' RA13T 'SCAL' PCG) '*' XVOLG  ;
  684.        dRA23T = ('EXTRAIRE' RA23T 'SCAL' PCD) '*' XVOLD  ;
  685.        gRA23T = ('EXTRAIRE' RA23T 'SCAL' PCG) '*' XVOLG  ;
  686.        HORROR = ('MAXIMUM' ('PROG' dRA13T gRA13T
  687.           dRA23T gRA23T) 'ABS') ;
  688.        ERRO = 1.0D-8 * fA12gd 'ABS' ;
  689.        LOGI2 = HORROR < ERRO ;
  690.        'SI' ('NON' LOGI2) ;
  691.           'MESSAGE' ('CHAINE' 'Erreur qdm_t ' HORROR ' > ' ERRO);
  692.        'FINSI' ;
  693.        LOGI1 = LOGI1 'ET' LOGI2;
  694. *
  695.        dRA14 = ('EXTRAIRE' RA14 'SCAL' PCD) '*' XVOLD ;
  696.        gRA14 = ('EXTRAIRE' RA14 'SCAL' PCG) '*' XVOLG ;
  697.        dRA24 = ('EXTRAIRE' RA24 'SCAL' PCD) '*' XVOLD ;
  698.        gRA24 = ('EXTRAIRE' RA24 'SCAL' PCG) '*' XVOLG ;
  699.        ERRO = 1.0D-7 * fA14gd 'ABS' ;
  700.        HORROR = ('MAXIMUM' ('PROG' gRA14 (dRA14 '-' fA14gd)
  701.           (gRA24 '+' fA24gd) dRA24) 'ABS') ;
  702.        LOGI2 = HORROR < ERRO ;
  703.        'SI' ('NON' LOGI2) ;
  704.           'MESSAGE' ('CHAINE' 'Erreur energy ' HORROR ' > ' ERRO);
  705.        'FINSI' ;
  706.        LOGI1 = LOGI1 'ET' LOGI2;
  707. *
  708.        'SI' ('NON' LOGI1);
  709.           'MESSAGE' ;
  710.           'MESSAGE' 'OPERATEUR KONV';
  711.           'MESSAGE' ('CHAINE' METO);
  712.           'MESSAGE' ;
  713.           'ERREUR'  ;
  714.        'FINSI' ;
  715.  
  716.  
  717.    'FIN' BLMETO ;       
  718.  
  719. ****************************************************
  720. ****************************************************
  721. ******** Fin boucle sur les angles         *********
  722. ****************************************************
  723. ****************************************************
  724.  
  725. 'FIN' BLOC;
  726.  
  727. 'FIN' ;
  728.  
  729.  
  730.  
  731.  
  732.  
  733.  
  734.  
  735.  
  736.  
  737.  
  738.  
  739.  
  740.  
  741.  
  742.  
  743.  
  744.  
  745.  

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