Download konv_resi_dem3D_stationaryshock_12.dgibi
* fichier : konv_resi_dem3D_stationaryshock_12.dgibi
***********************************************************
**** APPROCHE VF "Cell-Centred Formulation" pour la ****
**** solution des ****
**** Equations d'Euler pour un gaz parfait ****
**** Approche DEM pour la combustion ****
**** OPERATEURS PRET, KONV ****
**** ****
**** 3D ****
**** ****
**** Stationary shock (phase 1 '-' vs phase 2) ****
**** ****
**** Methodes: VLH, SS ****
**** ****
**** A. BECCANTINI DM2S/SFME/LTMF DECEMBRE 2010 ****
***********************************************************
'OPTION' 'DIME' 3 'ELEM' 'CUB8' 'ECHO' 0 'TRAC' 'X' ;
*
*** GRAPH
*
GRAPH = FAUX ;
* GRAPH = VRAI ;
**** Cas homogene : la table proprieté de gaz
*************************************************
**** The table for the properties of the gas ****
*************************************************
*
PGAS = 'TABLE' ;
*
**** Order of the polynomial order for cv = cv(T)
* For T > TMAX, cv(T) = cv(Tmax)
*
PGAS . 'TMAX' = 6000.0 ;
PGAS . 'NORD' = 4 ;
*
**** Species involved in the mixture (before or after
* the chemical reaction)
*
PGAS . 'SPECIES' = 'MOTS' 'H2 ' 'O2 ' 'H2O ' 'N2 ' ;
*
*
**** Coefficient of the chemical reaction.
* Note that for the first species this coefficient should be positive
* Normal, we take it equal to 1.
*
* H2 '+' 0.5 O2 ---> H2O
*
PGAS . 'CHEMCOEF' = 'PROG' 1.0 0.5 -1.0 0.0 ;
*
**** Mass fraction of the first species before and after the combustion
* Final mass fractions of the species with positive coefficients.
* Final mass fractions of the species with non-positive coefficient.
* The mass fraction of the last species is not given.
*
PGAS . 'MASSFRA' = 'PROG' 0.285219E-01 0.964039E-11 0.765104E-10
0.127442E-10 ;
*
**** Coef with the gas properties
*
PGAS . 'H2 ' = 'TABLE' ;
PGAS . 'H2O ' = 'TABLE' ;
PGAS . 'N2 ' = 'TABLE' ;
PGAS . 'O2 ' = 'TABLE' ;
*
**** Runiv (J/mole/K)
*
PGAS . 'RUNIV' = 8.31441 ;
*
**** W (kg/mole). Gas constant (J/kg/K = Runiv/W)
*
PGAS . 'H2 ' . 'W' = 2. * 1.00797E-3 ;
PGAS . 'O2 ' . 'W' = 2. * 15.9994E-3 ;
PGAS . 'H2O ' . 'W' = (PGAS . 'H2 ' . 'W' ) '+'
(0.5 * (PGAS . 'O2 ' . 'W' )) ;
PGAS . 'N2 ' . 'W' = 2 * 14.0067E-3 ;
*
**** Polynomial coefficients
*
PGAS . 'H2 ' . 'A' = 'PROG' 9834.91866 0.54273926 0.000862203836
-2.37281455E-07 1.84701105E-11 ;
PGAS . 'H2O ' . 'A' = 'PROG' 1155.95625 0.768331151 -5.73129958E-05
-1.82753232E-08 2.44485692E-12 ;
PGAS . 'N2 ' . 'A' = 'PROG' 652.940766 0.288239099 -7.80442298E-05
8.78233606E-09 -3.05514485E-13 ;
PGAS . 'O2 ' . 'A' = 'PROG' 575.012333 0.350522002 -0.000128294865
2.33636971E-08 -1.53304905E-12;
*
**** Formation enthalpies (energies) at 0K (J/Kg)
*
PGAS . 'H2 ' . 'H0K' = -4.195D6 ;
PGAS . 'H2O ' . 'H0K' = -1.395D7 ;
PGAS . 'N2 ' . 'H0K' = -2.953D5 ;
PGAS . 'O2 ' . 'H0K' = -2.634D5 ;
***************************
***** DOMAINE SPATIAL ****
***************************
A1 = 0.0D0 0.0D0 0.0;
A2 = 1.0D0 0.0D0 0.0;
A3 = 2.0D0 0.0D0 0.0;
A4 = 2.0D0 1.0D0 0.0;
A5 = 1.0D0 1.0D0 0.0;
A6 = 0.0D0 1.0D0 0.0;
L12 = A1 'DROIT' 1 A2;
L23 = A2 'DROIT' 1 A3;
L34 = A3 'DROIT' 1 A4;
L45 = A4 'DROIT' 1 A5;
L56 = A5 'DROIT' 1 A6;
L61 = A6 'DROIT' 1 A1;
L25 = A2 'DROIT' 1 A5;
BAS1 = 'DALL' L12 L25 L56 L61
'PLANE';
DOM1 = 'VOLUME' BAS1 'TRANSLATION' 1 (0.0 0.0 1.0) ;
BAS2 = 'DALL' L23 L34 L45 ('INVERSE' L25)
'PLANE';
DOM2 = 'VOLUME' BAS2 'TRANSLATION' 1 (0.0 0.0 1.0) ;
DOM10 = DOM1 ;
DOM20 = DOM2 ;
DOMTOT = DOM1 'ET' DOM2 ;
'ELIMINATION' 0.0001 DOMTOT ;
*
*** Point ou on controlle la consistence
*
P10 = 1.0 0.5 0.5 ;
$DOMTOT = 'MODELISER' DOMTOT 'EULER';
$DOM1 = 'MODELISER' DOM1 'EULER';
$DOM2 = 'MODELISER' DOM2 'EULER';
TDOMTOT = 'DOMA' $DOMTOT 'VF';
TDOM1 = 'DOMA' $DOM1 'VF';
TDOM2 = 'DOMA' $DOM2 'VF';
MDOM1 = TDOM1 . 'QUAF' ;
MDOM2 = TDOM2 . 'QUAF' ;
MDOMTOT = TDOMTOT . 'QUAF' ;
'ELIM' (MDOMTOT 'ET' MDOM1 'ET' MDOM2) 1.E-6 ;
'SI' GRAPH;
'TRACER' (('DOMA' $DOMTOT 'MAILLAGE') 'ET'
('DOMA' $DOMTOT 'FACEL') 'ET' P1) 'TITRE' 'Domaine et FACEL';
'FINSI' ;
*
EPS = 1.0D-16 ;
K0 = 45.2 ;
*
*******************************************
*** Etats gauche et droite ****************
*******************************************
*
* There are 3 fluxes to control :
* 1 -> 2 can be checked by considering an exact reactive shock
* 2 -> 1 can be checked by considering an exact reactive shock
* 1 -> 1 can be checked via consistency (or via an exact non-reactive
* or by comparing the results with a non-reactive one)
* 2 -> 2 can be checked via consistency (or via an exact non-reactive
* or by comparing the results with a non-reactive one)
*
* Case 1 -> 1
* Constant state
*
epsv = 0.0 ;
tg1 = 327.0 ;
td1 = 327.0 ;
pg1 = 139000.0 ;
pd1 = 139000.0 ;
ung1 = 45.2 '-' epsv ;
und1 = 45.2 '-' epsv ;
utg1 = 0.0 ;
utd1 = 0.0 ;
uvg1 = 0.0 ;
uvd1 = 0.0 ;
* utg1 = 0.0 ;
* utd1 = 0.0 ;
tg2 = 2518.582685186033 ;
td2 = 2518.582685186033 ;
pg2 = 125279.9833247705 ;
pd2 = 125279.9833247705 ;
ung2 = 329.121088555842 '-' epsv ;
und2 = 329.121088555842 '-' epsv ;
utg2 = 0.0 ;
utd2 = 0.0 ;
uvg2 = 0.0 ;
uvd2 = 0.0 ;
TN1 = 'PROG' tg1 td1 ;
TN2 = 'PROG' tg2 td2 ;
PN1 = 'PROG' pg1 pd1 ;
PN2 = 'PROG' pg2 pd2 ;
UN1 = 'PROG' ung1 und1 ;
UN2 = 'PROG' ung2 und2 ;
UT1 = 'PROG' utg1 utd1 ;
UT2 = 'PROG' utg2 utd2 ;
UV1 = 'PROG' uvg1 uvd1 ;
UV2 = 'PROG' uvg2 uvd2 ;
ALPHA1 = 'PROG' 1.0 0.0 ;
ALPHA2 = ('PROG' ('DIME' ALPHA1) '*' 1.0) '-' ALPHA1 ;
*
*
* ETHER = int_0^T cv(T') dT' T < TMAX
* = int_0^TMAX cv(T') dT' '+'
* cv(TMAX) T >= TMAX
*
ESP1 = 'EXTRAIRE' (PGAS . 'SPECIES') 1 ;
* DY1 = y_i - y_f for the species 1
DY1 = (('EXTRAIRE' (PGAS . 'MASSFRA') 1) '-'
('EXTRAIRE' (PGAS . 'MASSFRA') 2)) ;
COEF1 = ('EXTRAIRE' (PGAS . 'CHEMCOEF') 1) '*'
(PGAS . ESP1 . 'W') ;
YFINPH1 = 1.0 ;
YFINPH2 = 1.0 ;
'SI' (COEF1 > 0) ;
YPH2 = 'EXTRAIRE' (PGAS . 'MASSFRA') 2 ;
YPH1 = YPH2 '+' DY1 ;
'SINON' ;
YPH1 = 'EXTRAIRE' (PGAS . 'MASSFRA') 2 ;
YPH2 = YPH1 '-' DY1 ;
'FINSI' ;
YFINPH1 = YFINPH1 '-' YPH1 ;
YFINPH2 = YFINPH2 '-' YPH2 ;
PRYPH1 = 'PROG' YPH1 ;
PRYPH2 = 'PROG' YPH2 ;
'REPETER' BLESP (('DIME' (PGAS . 'SPECIES')) '-' 2) ;
ESP = 'EXTRAIRE' (PGAS . 'SPECIES') (&BLESP '+' 1) ;
COEF = ('EXTRAIRE' (PGAS . 'CHEMCOEF') (&BLESP '+' 1))
'*' (PGAS . ESP . 'W') ;
DY = (DY1 * (COEF '/' COEF1)) ;
'SI' (COEF > 0) ;
YPH2 = 'EXTRAIRE' (PGAS . 'MASSFRA') (&BLESP '+' 2) ;
YPH1 = YPH2 '+' DY ;
'SINON' ;
YPH1 = 'EXTRAIRE' (PGAS . 'MASSFRA') (&BLESP '+' 2) ;
YPH2 = YPH1 '-' DY ;
'FINSI' ;
PRYPH1 = PRYPH1 'ET' ('PROG' YPH1) ;
PRYPH2 = PRYPH2 'ET' ('PROG' YPH2) ;
YFINPH1 = YFINPH1 '-' YPH1 ;
YFINPH2 = YFINPH2 '-' YPH2 ;
'FIN' BLESP ;
PRYPH1 = PRYPH1 'ET' ('PROG' YFINPH1) ;
PRYPH2 = PRYPH2 'ET' ('PROG' YFINPH2) ;
'LISTE' PRYPH1 ;
'LISTE' PRYPH2 ;
*
TMAX = ('PROG' ('DIME' ALPHA1) '*' (PGAS . 'TMAX')) ;
* TCAL1 = MIN TN1, TMAX
TCAL1 = 0.5D0 '*' ((TMAX '+' TN1) '-' ('ABS' (TN1 '-' TMAX))) ;
DTN1 = TN1 '-' TCAL1 ;
* TCAL1 = MIN TN1, TMAX
TCAL2 = 0.5D0 '*' ((TMAX '+' TN2) '-' ('ABS' (TN2 '-' TMAX))) ;
DTN2 = TN2 '-' TCAL2 ;
*
* Internal energy (J/kg in SI)
*
ETHER1 = ('PROG' ('DIME' ALPHA1) '*' 0.0) ;
CV1 = ('PROG' ('DIME' ALPHA1) '*' 0.0) ;
ETHER2 = ('PROG' ('DIME' ALPHA1) '*' 0.0) ;
CV2 = ('PROG' ('DIME' ALPHA1) '*' 0.0) ;
FUNTN1 = ('PROG' ('DIME' ALPHA1) '*' 1.0) ;
FUNTN2 = ('PROG' ('DIME' ALPHA1) '*' 1.0) ;
'REPETER' BLPO ((PGAS . 'NORD') '+' 1) ;
'REPETER' BLESP ('DIME' (PGAS . 'SPECIES')) ;
ESP = 'EXTRAIRE' (PGAS . 'SPECIES') &BLESP ;
YCEL1 = 'EXTRAIRE' PRYPH1 &BLESP ;
YCEL2 = 'EXTRAIRE' PRYPH2 &BLESP ;
AA = 'EXTRAIRE' (PGAS . ESP . 'A') &BLPO ;
DCV1 = (AA * YCEL1 * FUNTN1) ;
DCV2 = (AA * YCEL2 * FUNTN2) ;
CV1 = CV1 '+' DCV1 ;
CV2 = CV2 '+' DCV2 ;
ETHER1 = ETHER1 '+' (DCV1 * TCAL1 '/' (&BLPO)) ;
ETHER2 = ETHER2 '+' (DCV2 * TCAL2 '/' (&BLPO)) ;
'FIN' BLESP ;
FUNTN1 = FUNTN1 '*' TCAL1 ;
FUNTN2 = FUNTN2 '*' TCAL2 ;
'FIN' BLPO ;
ETHER1 = ETHER1 '+' (CV1 '*' DTN1) ;
ETHER2 = ETHER2 '+' (CV2 '*' DTN2) ;
*
* Formation energy/enthalpy (J/kg in SI) and gas constant (J/kg/K)
*
EFORM1 = 0.0 ;
EFORM2 = 0.0 ;
RGAS1 = 0.0 ;
RGAS2 = 0.0 ;
'REPETER' BLESP ('DIME' (PGAS . 'SPECIES')) ;
ESP = 'EXTRAIRE' (PGAS . 'SPECIES') &BLESP ;
YCEL1 = 'EXTRAIRE' PRYPH1 &BLESP ;
YCEL2 = 'EXTRAIRE' PRYPH2 &BLESP ;
EFORM1 = EFORM1 '+' (YCEL1 * (PGAS . ESP . 'H0K')) ;
EFORM2 = EFORM2 '+' (YCEL2 * (PGAS . ESP . 'H0K')) ;
RGAS1 = RGAS1 '+' (YCEL1 * (PGAS . 'RUNIV') '/'
(PGAS . ESP . 'W')) ;
RGAS2 = RGAS2 '+' (YCEL2 * (PGAS . 'RUNIV') '/'
(PGAS . ESP . 'W')) ;
'FIN' BLESP ;
EFORM1 = 'PROG' ('DIME' ETHER1) '*' EFORM1 ;
EFORM2 = 'PROG' ('DIME' ETHER2) '*' EFORM2 ;
'LISTE' rgas1 ;
'LISTE' rgas2 ;
*
* Computation of the conservative variables
*
RN1 = PN1 '/' (RGAS1 '*' TN1) ;
RN2 = PN2 '/' (RGAS2 '*' TN2) ;
GNX1 = RN1 * UN1 ;
GNY1 = RN1 * UT1 ;
GNX2 = RN2 * UN2 ;
GNY2 = RN2 * UT2 ;
ECIN1 = 0.5D0 '*' ('PROG'
((ung1 * ung1) '+' (utg1 * utg1) '+' (uvg1 * uvg1))
((und1 * und1) '+' (utd1 * utd1) '+' (uvd1 * uvd1))
);
ECIN2 = 0.5D0 '*' ('PROG'
((ung2 * ung2) '+' (utg2 * utg2) '+' (uvg2 * uvg2))
((und2 * und2) '+' (utd2 * utd2) '+' (uvd2 * uvd2)));
RETN1 = RN1 '*' (ETHER1 '+' ECIN1 '+' EFORM1) ;
RETN2 = RN2 '*' (ETHER2 '+' ECIN2 '+' EFORM2) ;
*
* Computation of the resi contribution in (n,t)
*
* State 1
retg1 = MAXI RETN1 ;
rhog1 = 'MAXIMUM' RN1 ;
fA11gd = ung1 '*' rhog1 ;
fA12gd = (fA11gd '*' ung1) '+' pg1 ;
fA13gd = fA11gd '*' utg1 ;
fA14gd = ung1 '*' (retg1 '+' pg1);
*
* State 2
retg2 = MAXI RETN2 ;
rhog2 = 'MAXIMUM' RN2 ;
fA21gd = ung2 '*' rhog2 ;
fA22gd = (fA21gd '*' ung2) '+' pg2 ;
fA23gd = fA21gd '*' utg2 ;
fA24gd = ung2 '*' (retg2 '+' pg2);
*
****************************************************
****************************************************
******** Boucle sur les angles *********
****************************************************
****************************************************
DANGLE = 360 '/' 7.15;
ANGLE = 0.0 ;
'REPETER' BLOC 8;
*
*** Rotation
*
ANGLE = ANGLE '+' DANGLE;
ORIG = 0.0D0 0.0D0 0.0 ;
P1 = 1.0 1.0 1.0 ;
'MESSAGE' ;
'MESSAGE' (CHAIN 'Angle de rotation= ' ANGLE);
'MESSAGE' ;
DOM1 = DOM10 'TOURNER' ANGLE ORIG P1 ;
DOM2 = DOM20 'TOURNER' ANGLE ORIG P1 ;
P1FAC = P10 'TOURNER' ANGLE ORIG P1 ;
DOMTOT = DOM1 ET DOM2;
'ELIMINATION' DOMTOT 1D-6;
$DOMTOT = 'MODELISER' DOMTOT 'EULER';
$DOM1 = 'MODELISER' DOM1 'EULER';
$DOM2 = 'MODELISER' DOM2 'EULER';
TDOMTOT = 'DOMA' $DOMTOT 'VF';
TDOM1 = 'DOMA' $DOM1 'VF';
TDOM2 = 'DOMA' $DOM2 'VF';
MDOM1 = TDOM1 . 'QUAF' ;
MDOM2 = TDOM2 . 'QUAF' ;
MDOMTOT = TDOMTOT . 'QUAF' ;
'ELIM' (MDOMTOT 'ET' MDOM1 'ET' MDOM2) 1.E-6 ;
GRALP1 = 'NOMC' ('MOTS' 'UX' 'UY' 'UZ') (TDOMTOT . 'XXNORMAF')
('MOTS' 'P1DX' 'P1DY' 'P1DZ') ;
*
VINF = 'MANUEL' 'CHPO' (TDOMTOT . 'CENTRE') 1 'SCAL' 100. ;
*
**** Redefinition de P1FAC dans $DOMTOT 'FACE'
*
P1FAC = ('DOMA' $DOMTOT 'FACE') 'POIN' 'PROC' P1FAC ;
XSURF = 'EXTRAIRE' ('DOMA' $DOMTOT 'XXSURFAC') P1FAC 'SCAL' ;
NX = 'EXTRAIRE' (TDOMTOT . 'MATROT') P1FAC 'MX' ;
NY = 'EXTRAIRE' (TDOMTOT . 'MATROT') P1FAC 'MY' ;
NZ = 'EXTRAIRE' (TDOMTOT . 'MATROT') P1FAC 'MZ' ;
TX = 'EXTRAIRE' (TDOMTOT . 'MATROT') P1FAC 'UX' ;
TY = 'EXTRAIRE' (TDOMTOT . 'MATROT') P1FAC 'UY' ;
TZ = 'EXTRAIRE' (TDOMTOT . 'MATROT') P1FAC 'UZ' ;
VX = 'EXTRAIRE' (TDOMTOT . 'MATROT') P1FAC 'RX' ;
VY = 'EXTRAIRE' (TDOMTOT . 'MATROT') P1FAC 'RY' ;
VZ = 'EXTRAIRE' (TDOMTOT . 'MATROT') P1FAC 'RZ' ;
GEOPC = ('DOMA' $DOMTOT 'FACEL') 'ELEM' 'APPUYE' 'LARGEMENT'
P1FAC ;
AA = 'CHANGER' GEOPC 'POI1' ;
PC1 = 'POIN' 1 AA ;
P2 = 'POIN' 1 ('DOMA' $DOMTOT 'CENTRE') ;
PC2 = 'POIN' 3 AA ;
* Si P3 n'existe pas, probleme en FACEL.
'SI' ('NEG' PC1 P2) ;
ORIENT = -1 ;
PCD = PC1 ;
PCG = PC2 ;
XVOLD = 'MAXIMUM' (TDOM1 . 'XXVOLUM') ;
XVOLG = 'MAXIMUM' (TDOM2 . 'XXVOLUM') ;
'SINON' ;
ORIENT = 1 ;
PCD = PC2 ;
PCG = PC1 ;
XVOLG = 'MAXIMUM' (TDOM1 . 'XXVOLUM') ;
XVOLD = 'MAXIMUM' (TDOM2 . 'XXVOLUM') ;
'FINSI' ;
uxg1 = ((ung1 '*' NX) '+' (utg1 '*' TX) '+' (uvg1 '*' VX))
'*' ORIENT ;
uyg1 = ((ung1 '*' NY) '+' (utg1 '*' TY) '+' (uvg1 '*' VY))
'*' ORIENT ;
uzg1 = ((ung1 '*' NZ) '+' (utg1 '*' TZ) '+' (uvg1 '*' VZ))
'*' ORIENT ;
uxg2 = ((ung2 '*' NX) '+' (utg2 '*' TX) '+' (uvg2 '*' VX))
'*' ORIENT ;
uyg2 = ((ung2 '*' NY) '+' (utg2 '*' TY) '+' (uvg2 '*' VY))
'*' ORIENT ;
uzg2 = ((ung2 '*' NZ) '+' (utg2 '*' TZ) '+' (uvg2 '*' VZ))
'*' ORIENT ;
uxd1 = ((und1 '*' NX) '+' (utd1 '*' TX) '+' (uvd1 '*' VX))
'*' ORIENT ;
uyd1 = ((und1 '*' NY) '+' (utd1 '*' TY) '+' (uvd1 '*' VY))
'*' ORIENT ;
uzd1 = ((und1 '*' NZ) '+' (utd1 '*' TZ) '+' (uvd1 '*' VZ))
'*' ORIENT ;
uxd2 = ((und2 '*' NX) '+' (utd2 '*' TX) '+' (uvd2 '*' VX))
'*' ORIENT ;
uyd2 = ((und2 '*' NY) '+' (utd2 '*' TY) '+' (uvd2 '*' VY))
'*' ORIENT ;
uzd2 = ((und2 '*' NZ) '+' (utd2 '*' TZ) '+' (uvd2 '*' VZ))
'*' ORIENT ;
* 'OPTION' DONN 5 ;
'SI' GRAPH;
'TRACER' (('DOMA' $DOMTOT 'MAILLAGE') 'ET'
('COULEUR' ('DOMA' $DOMTOT 'FACEL') 'ROUGE') 'ET' P1FAC)
'TITRE' 'Domaine et FACEL';
'FINSI' ;
MAILIM = 'DIFF' ('DOMA' $DOMTOT 'FACE') ('MANUEL' 'POI1' P1FAC) ;
***********************
**** Les CHPOINTs ****
***********************
CHRN1 = ('MANUEL' 'CHPO' (TDOMTOT . 'CENTRE') 1 'SCAL' RN1) ;
CHRN2 = ('MANUEL' 'CHPO' (TDOMTOT . 'CENTRE') 1 'SCAL' RN2) ;
CHVN1 = ('MANUEL' 'CHPO' (TDOM1 . 'CENTRE') 3 'UX' uxg1
'UY' uyg1 'UZ' uzg1) '+'
('MANUEL' 'CHPO' (TDOM2 . 'CENTRE') 3 'UX' uxd1
'UY' uyd1 'UZ' uzd1) ;
CHVN2 = ('MANUEL' 'CHPO' (TDOM1 . 'CENTRE') 3 'UX' uxg2
'UY' uyg2 'UZ' uzg2) '+'
('MANUEL' 'CHPO' (TDOM2 . 'CENTRE') 3 'UX' uxd2
'UY' uyd2 'UZ' uzd2) ;
CHGN1 = CHRN1 '*' CHVN1 ;
CHGN2 = CHRN2 '*' CHVN2 ;
CHRET1 = ('MANUEL' 'CHPO' (TDOMTOT . 'CENTRE') 1 'SCAL' RETN1) ;
CHRET2 = ('MANUEL' 'CHPO' (TDOMTOT . 'CENTRE') 1 'SCAL' RETN2) ;
CHTN1 = ('MANUEL' 'CHPO' (TDOMTOT . 'CENTRE') 1 'SCAL' TN1) ;
CHTN2 = ('MANUEL' 'CHPO' (TDOMTOT . 'CENTRE') 1 'SCAL' TN2) ;
CHAL1 = ('MANUEL' 'CHPO' (TDOMTOT . 'CENTRE') 1 'SCAL' ALPHA1) ;
CHAL2 = ('MANUEL' 'CHPO' (TDOMTOT . 'CENTRE') 1 'SCAL' ALPHA2) ;
R1 R2 V1 V2 P1 P2 T1 T2 = 'PRIM' 'DEM' PGAS CHAL1 CHAL2
(CHAL1 * CHRN1) (CHAL2 * CHRN2) (CHAL1 * CHGN1)
(CHAL2 * CHGN2) (CHAL1 * CHRET1) (CHAL2 * CHRET2)
CHTN1 CHTN2 EPS ;
GRADAL1 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 3
'P1DX' 0.0 'P1DY' 0.0 'P1DZ' 0.0 ;
GRADAL2 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 3
'P1DX' 0.0 'P1DY' 0.0 'P1DZ' 0.0 ;
LIMAL1 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 1
'P1' 0.0 ;
LIMAL2 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 1
'P1' 0.0 ;
GRADR1 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 3
'P1DX' 0.0 'P1DY' 0.0 'P1DZ' 0.0;
GRADR2 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 3
'P1DX' 0.0 'P1DY' 0.0 'P1DZ' 0.0;
LIMR1 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 1
'P1' 0.0 ;
LIMR2 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 1
'P1' 0.0 ;
GRADV1 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 9
'P1DX' 0.0 'P1DY' 0.0 'P1DZ' 0.0
'P2DX' 0.0 'P2DY' 0.0 'P2DZ' 0.0
'P3DX' 0.0 'P3DY' 0.0 'P3DZ' 0.0 ;
GRADV2 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 9
'P1DX' 0.0 'P1DY' 0.0 'P1DZ' 0.0
'P2DX' 0.0 'P2DY' 0.0 'P2DZ' 0.0
'P3DX' 0.0 'P3DY' 0.0 'P3DZ' 0.0 ;
LIMV1 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 3
'P1' 0.0 'P2' 0.0 'P3' 0.0 ;
LIMV2 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 3
'P1' 0.0 'P2' 0.0 'P3' 0.0 ;
GRADP1 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 3
'P1DX' 0.0 'P1DY' 0.0 'P1DZ' 0.0;
GRADP2 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 3
'P1DX' 0.0 'P1DY' 0.0 'P1DZ' 0.0;
LIMP1 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 1
'P1' 0.0 ;
LIMP2 = 'MANU' 'CHPO' (TDOMTOT . 'CENTRE') 1
'P1' 0.0 ;
*
**** L'operateur 'PRET'
*
CHFAL1 CHFAL2 CHFRN1 CHFRN2 CHFVN1 CHFVN2 CHFPN1 CHFPN2 =
'PRET' 'DEM' $DOMTOT
CHAL1 (0.0 * GRADAL1) LIMAL1
CHAL2 (0.0 * GRADAL2) LIMAL2
R1 (0.0 * GRADR1) LIMR1
R2 (0.0 * GRADR2) LIMR2
V1 (0.0 * GRADV1) LIMV1
V2 (0.0 * GRADV2) LIMV2
P1 (0.0 * GRADP1) LIMP1
P2 (0.0 * GRADP2) LIMP2 ;
*
SI FAUX ;
CHFRO1 = 'REDU' CHFRN1 geopc ;
CHFV1 = 'REDU' CHFVN1 geopc ;
CHFP1 = 'REDU' CHFPN1 geopc ;
CHFRO2 = 'REDU' CHFRN2 geopc ;
CHFV2 = 'REDU' CHFVN2 geopc ;
CHFP2 = 'REDU' CHFPN2 geopc ;
'LISTE' CHFRO1 ;
'LISTE' CHFV1 ;
'LISTE' CHFP1 ;
'LISTE' CHFRO2 ;
'LISTE' CHFV2 ;
'LISTE' CHFP2 ;
'OPTION' donn 5 ;
'FINSI' ;
*
**** L'operateur 'KONV'
*
'REPETER' BLMETO 3 ;
'SI' ('EGA' &BLMETO 1) ;
METO = 'SS' ;
'FINSI' ;
'SI' ('EGA' &BLMETO 2) ;
METO = 'VLH' ;
'FINSI' ;
'SI' ('EGA' &BLMETO 3) ;
METO = 'AUSMPUP' ;
'FINSI' ;
'MESSAGE' ('CHAINE' 'METO = ' METO) ;
LISTINC1 = ('MOTS' 'ALF1' 'RN1' 'RNX1' 'RUY1' 'RUZ1' 'RET1'
'ALF2' 'RN2' 'RUX2' 'RUY2' 'RUZ2' 'RET2') ;
SI ('EGA' METO 'AUSMPUP') ;
CHPRES DT SURF = 'KONV' 'VF' 'DEM' 'RESI' METO 'CONS'
$DOMTOT PGAS LISTINC1 CHAL1 CHAL2 CHFAL1 CHFAL2 CHFRN1 CHFRN2
CHFVN1 CHFVN2 CHFPN1 CHFPN2 K0 GRALP1 EPS MAILIM VINF VINF ;
'SINON' ;
CHPRES DT SURF = 'KONV' 'VF' 'DEM' 'RESI' METO 'CONS'
$DOMTOT PGAS LISTINC1 CHAL1 CHAL2 CHFAL1 CHFAL2 CHFRN1 CHFRN2
CHFVN1 CHFVN2 CHFPN1 CHFPN2 K0 GRALP1 EPS MAILIM ;
'FINSI' ;
* 'LISTE' ('CHAINE' 'FORMAT' '(E16.10)' ('EXTRAIRE' ETHER1 1)) ;
* 'LISTE' ('CHAINE' 'FORMAT' '(E16.10)' ('EXTRAIRE' ECIN1 1)) ;
* 'LISTE' ('CHAINE' 'FORMAT' '(E16.10)' ('EXTRAIRE' EFORM1 1)) ;
* 'LISTE' ('CHAINE' 'FORMAT' '(E16.10)' (pg1)) ;
* 'LISTE' ('CHAINE' 'FORMAT' '(E16.10)' ('EXTRAIRE' RN1 1)) ;
* 'LISTE' ('CHAINE' 'FORMAT' '(E16.10)' (ung1)) ;
* 'OPTION' DONN 5 ;
*
* We verify that the residum contribution of P1FAC is the one
* of a reactive stationary shock.
*
RA11 = 'EXCO' ('EXTRAIRE' LISTINC1 1) CHPRES ;
RA12 = 'EXCO' ('EXTRAIRE' LISTINC1 2) CHPRES ;
RA13X = 'EXCO' ('EXTRAIRE' LISTINC1 3) CHPRES ;
RA13Y = 'EXCO' ('EXTRAIRE' LISTINC1 4) CHPRES ;
RA13Z = 'EXCO' ('EXTRAIRE' LISTINC1 5) CHPRES ;
RA13N = (RA13X '*' NX) '+' (RA13Y * NY) '+' (RA13Z * NZ) ;
RA13T = (RA13X '*' TX) '+' (RA13Y * TY) '+' (RA13Z * TZ) ;
RA13V = (RA13X '*' VX) '+' (RA13Y * VY) '+' (RA13Z * VZ) ;
RA14 = 'EXCO' ('EXTRAIRE' LISTINC1 6) CHPRES ;
*
RA21 = 'EXCO' ('EXTRAIRE' LISTINC1 7) CHPRES ;
RA22 = 'EXCO' ('EXTRAIRE' LISTINC1 8) CHPRES ;
RA23X = 'EXCO' ('EXTRAIRE' LISTINC1 9) CHPRES ;
RA23Y = 'EXCO' ('EXTRAIRE' LISTINC1 10) CHPRES ;
RA23Z = 'EXCO' ('EXTRAIRE' LISTINC1 11) CHPRES ;
RA23N = (RA23X '*' NX) '+' (RA23Y * NY) '+' (RA23Z * NZ) ;
RA23T = (RA23X '*' TX) '+' (RA23Y * TY) '+' (RA23Z * TZ) ;
RA23V = (RA23X '*' VX) '+' (RA23Y * VY) '+' (RA23Z * VZ) ;
RA24 = 'EXCO' ('EXTRAIRE' LISTINC1 12) CHPRES ;
*
ERRO = 1D-7 ;
HORROR = 'MAXIMUM' RA11 'ABS' ;
LOGI1 = HORROR < ERRO ;
'SI' ('NON' LOGI1) ;
'MESSAGE' ('CHAINE' 'Erreur alpha ' HORROR ' > ' ERRO);
'FINSI' ;
HORROR = 'MAXIMUM' RA21 'ABS' ;
LOGI2 = HORROR < ERRO ;
'SI' ('NON' LOGI2) ;
'MESSAGE' ('CHAINE' 'Erreur alpha ' HORROR ' > ' ERRO);
'FINSI' ;
LOGI1 = LOGI1 'ET' LOGI2;
*
dRA12 = ('EXTRAIRE' RA12 'SCAL' PCD) '*' XVOLD ;
gRA12 = ('EXTRAIRE' RA12 'SCAL' PCG) '*' XVOLG ;
dRA22 = ('EXTRAIRE' RA22 'SCAL' PCD) '*' XVOLD ;
gRA22 = ('EXTRAIRE' RA22 'SCAL' PCG) '*' XVOLG ;
ERRO = 1.0D-8 * fA11gd 'ABS' ;
HORROR = ('MAXIMUM' ('PROG' dRA12 (gRA12 '+' fA11gd)
(dRA22 '-' fA21gd) gRA22) 'ABS') ;
LOGI2 = HORROR < ERRO ;
'SI' ('NON' LOGI2) ;
'MESSAGE' ('CHAINE' 'Erreur density ' HORROR ' > ' ERRO);
'FINSI' ;
LOGI1 = LOGI1 'ET' LOGI2;
*
dRA13N = ('EXTRAIRE' RA13N 'SCAL' PCD) '*' XVOLD ;
gRA13N = ('EXTRAIRE' RA13N 'SCAL' PCG) '*' XVOLG ;
dRA23N = ('EXTRAIRE' RA23N 'SCAL' PCD) '*' XVOLD ;
gRA23N = ('EXTRAIRE' RA23N 'SCAL' PCG) '*' XVOLG ;
ERRO = 1.0D-8 * fA12gd 'ABS' ;
HORROR = ('MAXIMUM' ('PROG' dRA13N (gRA13N '+' fA12gd)
(dRA23N '-' fA22gd) gRA23N) 'ABS') ;
LOGI2 = HORROR < ERRO ;
'SI' ('NON' LOGI2) ;
'MESSAGE' ('CHAINE' 'Erreur qdm_n ' HORROR ' > ' ERRO);
'FINSI' ;
LOGI1 = LOGI1 'ET' LOGI2;
*
dRA13T = ('EXTRAIRE' RA13T 'SCAL' PCD) '*' XVOLD ;
gRA13T = ('EXTRAIRE' RA13T 'SCAL' PCG) '*' XVOLG ;
dRA23T = ('EXTRAIRE' RA23T 'SCAL' PCD) '*' XVOLD ;
gRA23T = ('EXTRAIRE' RA23T 'SCAL' PCG) '*' XVOLG ;
HORROR = ('MAXIMUM' ('PROG' dRA13T gRA13T
dRA23T gRA23T) 'ABS') ;
ERRO = 1.0D-8 * fA12gd 'ABS' ;
LOGI2 = HORROR < ERRO ;
'SI' ('NON' LOGI2) ;
'MESSAGE' ('CHAINE' 'Erreur qdm_t ' HORROR ' > ' ERRO);
'FINSI' ;
LOGI1 = LOGI1 'ET' LOGI2;
*
dRA14 = ('EXTRAIRE' RA14 'SCAL' PCD) '*' XVOLD ;
gRA14 = ('EXTRAIRE' RA14 'SCAL' PCG) '*' XVOLG ;
dRA24 = ('EXTRAIRE' RA24 'SCAL' PCD) '*' XVOLD ;
gRA24 = ('EXTRAIRE' RA24 'SCAL' PCG) '*' XVOLG ;
ERRO = 1.0D-7 * fA14gd 'ABS' ;
HORROR = ('MAXIMUM' ('PROG' dRA14 (gRA14 '+' fA14gd)
(dRA24 '-' fA24gd) gRA24) 'ABS') ;
LOGI2 = HORROR < ERRO ;
'SI' ('NON' LOGI2) ;
'MESSAGE' ('CHAINE' 'Erreur energy ' HORROR ' > ' ERRO);
'FINSI' ;
LOGI1 = LOGI1 'ET' LOGI2;
*
'SI' ('NON' LOGI1);
'MESSAGE' ;
'MESSAGE' 'OPERATEUR KONV';
'MESSAGE' ('CHAINE' METO);
'MESSAGE' ;
'ERREUR' ;
'FINSI' ;
'FIN' BLMETO ;
****************************************************
****************************************************
******** Fin boucle sur les angles *********
****************************************************
****************************************************
'FIN' BLOC;
'FIN' ;
© Cast3M 2003 - All rights reserved.
Disclaimer