Download konv_ther_cons2.dgibi
* fichier : konv_ther_cons2.dgibi
***********************************************************
**** APPROCHE VF "Cell-Centred Formulation" pour la ****
**** solution des ****
**** Equations d'Euler pour un gaz parfait ****
**** OPERATEURS PRET, KONV ****
**** ****
**** Cas gaz multi-especes "thermally perfect" ****
**** Consistence ****
**** Splitting des scalaires passifs ****
**** ****
**** Methodes: VLH, Colella ****
**** ****
**** A. BECCANTINI DRN/DMT/SEMT/LTMF FEVRIER 2000 ****
***********************************************************
'OPTION' 'DIME' 2 ;
'OPTION' 'ELEM' QUA4 ;
'OPTION' 'ECHO' 0 ;
'OPTION' 'TRAC' 'X' ;
*
*** GRAPH
*
GRAPH = FAUX ;
* GRAPH = VRAI ;
***********************
**** LA TABLE PGAZ ****
***********************
PGAZ = 'TABLE' ;
PGAZ2 = 'TABLE' ;
*
**** Ordre des polynoms cv_i
*
PGAZ . 'NORD' = 1 ;
PGAZ2 . 'NORD' = 1 ;
*
**** Especes qui sont dans les equations d'Euler
*
PGAZ . 'ESPEULE' = 'MOTS' 'H2 ' 'O2 ' 'H2O ' ;
PGAZ2 . 'ESPEULE' = 'MOTS' 'H2 ' 'O2 ' 'H2O ' ;
*
**** Espece qui n'y est pas
*
PGAZ . 'ESPNEULE' = 'N2 ';
PGAZ2 . 'ESPNEULE' = 'N2 ';
*
PGAZ . 'H2 ' = 'TABLE' ;
PGAZ . 'H2O ' = 'TABLE' ;
PGAZ . 'N2 ' = 'TABLE' ;
PGAZ . 'O2 ' = 'TABLE' ;
PGAZ2 . 'H2 ' = 'TABLE' ;
PGAZ2 . 'H2O ' = 'TABLE' ;
PGAZ2 . 'N2 ' = 'TABLE' ;
PGAZ2 . 'O2 ' = 'TABLE' ;
*
**** R (J/Kg/K)
*
PGAZ . 'H2 ' . 'R' = 4130.0 ;
PGAZ . 'H2O ' . 'R' = 461.4 ;
PGAZ . 'N2 ' . 'R' = 296.8 ;
PGAZ . 'O2 ' . 'R' = 259.8 ;
PGAZ2 . 'H2 ' . 'R' = 4130.0 ;
PGAZ2 . 'H2O ' . 'R' = 461.4 ;
PGAZ2 . 'N2 ' . 'R' = 296.8 ;
PGAZ2 . 'O2 ' . 'R' = 259.8 ;
*
**** Regressions polynomials
*
PGAZ . 'H2 ' . 'A' = 'PROG' 9834.91866 0.54273926 ;
PGAZ . 'H2O ' . 'A' = 'PROG' 1155.95625 0.768331151 ;
PGAZ . 'N2 ' . 'A' = 'PROG' 652.940766 0.288239099 ;
PGAZ . 'O2 ' . 'A' = 'PROG' 575.012333 0.350522002 ;
PGAZ2 . 'H2 ' . 'A' = 'PROG' 9834.91866 0.54273926 ;
PGAZ2 . 'H2O ' . 'A' = 'PROG' 1155.95625 0.768331151 ;
PGAZ2 . 'N2 ' . 'A' = 'PROG' 652.940766 0.288239099 ;
PGAZ2 . 'O2 ' . 'A' = 'PROG' 575.012333 0.350522002 ;
*
**** "Enthalpies" (ou energies) de formations a OK (J/Kg)
* Note: ce sont des entites numeriques
PGAZ . 'H2 ' . 'H0K' = -4.195D6 ;
PGAZ . 'H2O ' . 'H0K' = -1.395D7 ;
PGAZ . 'N2 ' . 'H0K' = -2.953D5 ;
PGAZ . 'O2 ' . 'H0K' = -2.634D5 ;
PGAZ2 . 'H2 ' . 'H0K' = -4.195D6 ;
PGAZ2 . 'H2O ' . 'H0K' = -1.395D7 ;
PGAZ2 . 'N2 ' . 'H0K' = -2.953D5 ;
PGAZ2 . 'O2 ' . 'H0K' = -2.634D5 ;
* PGAZ2: transport de scalaires passifs
PGAZ2 . 'SCALPASS' = 'MOTS' 'TUTU' 'TATA' ;
***************************
***** DOMAINE SPATIAL ****
***************************
A1 = 0.0D0 0.0D0;
A2 = 1.0D0 0.0D0;
A3 = 2.0D0 0.0D0;
A4 = 2.0D0 1.0D0;
A5 = 1.0D0 1.0D0;
A6 = 0.0D0 1.0D0;
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;
DOM10 = 'DALL' L12 L25 L56 L61
'PLANE';
DOM20 = 'DALL' L23 L34 L45 ('INVERSE' L25)
'PLANE';
*
*** Point ou on controlle la consistence
*
P10 = 1.0 0.5;
*
*** Etats gauche et droite
*
rogd = 1.412 ;
pgd = 1.17D6 ;
ungd = 8.34 ;
utgd = 18.2 ;
yh2gd = 0.1 ;
yh2Ogd = 0.2 ;
yN2gd = 0.15 ;
yO2gd = 1 '-' (yh2gd '+' yh2Ogd '+' yn2gd) ;
tutugd = 14.9 ;
tatagd = 23.8 ;
*
**** tg
*
tgd = pgd '/' (rogd *
((yh2gd * (PGAZ . 'H2 ' . 'R')) '+'
(yh2Ogd * (PGAZ . 'H2O ' . 'R')) '+'
(yn2gd * (PGAZ . 'N2 ' . 'R')) '+'
(yo2gd * (PGAZ . 'O2 ' . 'R')))) ;
*
**** gamg, gamd , htg , htd
*
*
A0H2 = 'EXTRAIRE' (PGAZ . 'H2 ' . 'A') 1 ;
A1H2 = 'EXTRAIRE' (PGAZ . 'H2 ' . 'A') 2 ;
A0H2O= 'EXTRAIRE' (PGAZ . 'H2O ' . 'A') 1 ;
A1H2O= 'EXTRAIRE' (PGAZ . 'H2O ' . 'A') 2 ;
A0N2 = 'EXTRAIRE' (PGAZ . 'N2 ' . 'A') 1 ;
A1N2 = 'EXTRAIRE' (PGAZ . 'N2 ' . 'A') 2 ;
A0O2 = 'EXTRAIRE' (PGAZ . 'O2 ' . 'A') 1 ;
A1O2 = 'EXTRAIRE' (PGAZ . 'O2 ' . 'A') 2 ;
T = tgd ;
T2 = T * T ;
ETGD = (yh2gd *
((A0H2 * T) '+' (A1H2 * T2 '/' 2)));
ETGD = ETGD '+' (yH2Ogd *
((A0H2O * T) '+' (A1H2O * T2 '/' 2)));
ETGD = ETGD '+' (yN2gd *
((A0N2 * T) '+' (A1N2 * T2 '/' 2)));
ETGD = ETGD '+' (yO2gd *
((A0O2 * T) '+' (A1O2 * T2 '/' 2))) ;
ecingd = 0.5D0 '*' ((ungd '*' ungd) '+' (utgd '*' utgd));
retgd = rogd '*' (ETGd '+' ecingd) ;
*
*** flux en (n,t)
*
f1gd = ungd '*' rogd ;
f2gd = (rogd '*' (ungd '*' ungd)) '+' pgd ;
f3gd = rogd '*' ungd '*' utgd ;
f4gd = ungd '*' (retgd '+' pgd);
f5gd = ungd '*' rogd '*' yh2gd ;
f6gd = ungd '*' rogd '*' yo2gd ;
f7gd = ungd '*' rogd '*' yh2ogd ;
ftutugd = ungd '*' rogd '*' tutugd ;
ftatagd = ungd '*' rogd '*' tatagd ;
****************************************************
****************************************************
******** Boucle sur les angles *********
****************************************************
****************************************************
DANGLE = 360 '/' 7.15;
ANGLE = 11.3 ;
'REPETER' BLMET 2 ;
SI (&BLMET 'EGA' 1) ;
METO = 'VLH' ;
'FINSI' ;
SI (&BLMET 'EGA' 2) ;
METO = 'SS' ;
'FINSI' ;
'REPETER' BLOC 8;
*
*** Rotation
*
ANGLE = ANGLE '+' DANGLE;
ORIG = 0.0D0 0.0D0;
MESSAGE;
MESSAGE (CHAIN 'Angle de rotation= ' ANGLE);
MESSAGE;
DOM1 = DOM10 'TOURNER' ANGLE ORIG;
DOM2 = DOM20 'TOURNER' ANGLE ORIG;
P1 = P10 'TOURNER' ANGLE ORIG;
'ELIMINATION' (DOM1 ET DOM2) 1D-6;
DOMTOT = DOM1 ET DOM2;
$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' ;
*
**** Redefinition de P1 dans $DOMTOT 'FACE'
*
P1 = ('DOMA' $DOMTOT 'FACE') 'POIN' 'PROC' P1;
***********************
**** Les CHPOINTs ****
***********************
uxgd = (ungd '*' ('COS' ANGLE)) '-' (utgd '*' ('SIN' ANGLE));
uygd = (ungd '*' ('SIN' ANGLE)) '+' (utgd '*' ('COS' ANGLE));
RN1 = 'KCHT' $DOM1 'SCAL' 'CENTRE' rogd;
RN2 = 'KCHT' $DOM2 'SCAL' 'CENTRE' rogd;;
RN = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' (RN1 'ET' RN2);
VN1 = 'KCHT' $DOM1 'VECT' 'CENTRE' (uxgd uygd);
VN2 = 'KCHT' $DOM2 'VECT' 'CENTRE' (uxgd uygd);
VN = 'KCHT' $DOMTOT 'VECT' 'CENTRE' (VN1 'ET' VN2);
PN1 = 'KCHT' $DOM1 'SCAL' 'CENTRE' pgd;
PN2 = 'KCHT' $DOM2 'SCAL' 'CENTRE' pgd;
PN = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' (PN1 'ET' PN2);
* attention: PGAZ . 'ESPEULE' = 'MOTS' 'H2 ' 'O2 ' 'H2O ' ;
YH2 = ('KCHT' $DOM1 'SCAL' 'CENTRE' 'COMP' 'H2 ' yh2gd) 'ET'
('KCHT' $DOM2 'SCAL' 'CENTRE' 'COMP' 'H2 ' yh2gd) ;
YO2 = ('KCHT' $DOM1 'SCAL' 'CENTRE' 'COMP' 'O2 ' yo2gd) 'ET'
('KCHT' $DOM2 'SCAL' 'CENTRE' 'COMP' 'O2 ' yo2gd) ;
YH2O = ('KCHT' $DOM1 'SCAL' 'CENTRE' 'COMP' 'H2O ' yh2ogd) 'ET'
('KCHT' $DOM2 'SCAL' 'CENTRE' 'COMP' 'H2O ' yh2ogd) ;
YN = YH2 'ET' YO2 'ET' YH2O ;
*
SCAN1 = 'KCHT' $DOM1 'VECT' 'CENTRE'
'COMP' (PGAZ2 . 'SCALPASS') (tutugd tatagd) ;
SCAN2 = 'KCHT' $DOM2 'VECT' 'CENTRE'
'COMP' (PGAZ2 . 'SCALPASS') (tutugd tatagd) ;
SCAN = 'KCHT' $DOMTOT 'VECT' 'CENTRE'
'COMP' (PGAZ2 . 'SCALPASS') (SCAN1 'ET' SCAN2);
***************************
**** L'operateur PRET****
***************************
ORDESP = 1;
ORDTEM = 1;
ROF VITF PF YF = 'PRET' 'PERFTEMP' ORDESP ORDTEM
$DOMTOT PGAZ RN VN PN YN ;
ROF2 VITF2 PF2 YF2 SCAF = 'PRET' 'PERFTEMP' ORDESP ORDTEM
$DOMTOT PGAZ2 RN VN PN YN SCAN ;
GEOP1 = ('DOMA' $DOMTOT 'FACEL') 'ELEM' 'APPUYE' 'LARGEMENT' P1;
GEOP2 = ('DOMA' $DOMTOT 'FACE') 'ELEM' 'APPUYE' 'LARGEMENT' P1;
REFGEOP1 = 'REDU' VITF GEOP2;
*
**** Orientation de la normal n de castem par raport a la
* notre; t est par consequence
*
NCOS = 'EXTRAIRE' REFGEOP1 'NX' 1 1 1;
NSIN = 'EXTRAIRE' REFGEOP1 'NY' 1 1 1;
'SI' (('ABS' NCOS) > ('ABS' NSIN));
ORIENT = ('COS' ANGLE) '/' NCOS;
'SINON';
ORIENT = ('SIN' ANGLE) '/' NSIN;
'FINSI' ;
ORIENT = 'SIGNE' 'FLOTTANT' ORIENT;
**************************
**************************
**** 'KONV', TEST VLH ****
**************************
**************************
LISTINCO = 'MOTS' 'RN' 'RUX' 'RUY' 'RETN' 'H2' 'O2' 'H2O' ;
TFLUX DELTAT = 'KONV' 'VF' 'PERFTEMP' 'FLUX' METO
$DOMTOT PGAZ LISTINCO ROF VITF PF YF ;
*
**** Les flux aux interfaces sont dans de
* CHPOINT FACE
*
FLUX1 = 'EXCO' 'RN' TFLUX ;
FLUX2X = 'EXCO' 'RUX' TFLUX ;
FLUX2Y = 'EXCO' 'RUY' TFLUX ;
FLUX3 = 'EXCO' 'RETN' TFLUX ;
FLUX4 = 'EXCO' ('MOTS' 'H2' 'O2' 'H2O') TFLUX
('MOTS' 'H2' 'O2' 'H2O') ;
FLUX2N = (FLUX2X '*' ('COS' ANGLE)) '+' (FLUX2Y * ('SIN' ANGLE));
FLUX2T = (FLUX2Y '*' ('COS' ANGLE)) '-' (FLUX2X * ('SIN' ANGLE));
*
**** Test 'KONV'
*
f1 = 'EXTRAIRE' FLUX1 'SCAL' P1;
ERRO = 'ABS' (1D-12 '*' f1gd);
LOGI1 = ('ABS' ((f1 '*' ORIENT) '-' f1gd)) < ERRO;
f2 = 'EXTRAIRE' FLUX2N 'SCAL' P1;
ERRO = 'ABS' (1D-12 '*' f2gd);
LOGI2 = ('ABS' ((f2 '*' ORIENT) '-' f2gd)) < ERRO;
LOGI1 = LOGI1 'ET' LOGI2;
f3 = 'EXTRAIRE' FLUX2T 'SCAL' P1;
ERRO = 'ABS' (1D-11 '*' f3gd) ;
LOGI2 = ('ABS' ((f3 '*' ORIENT) '-' f3gd)) < ERRO;
LOGI1 = LOGI1 'ET' LOGI2;
f4 = 'EXTRAIRE' FLUX3 'SCAL' P1;
ERRO = 'ABS' (1D-12 '*' f4gd) ;
LOGI2 = ('ABS' ((f4 '*' ORIENT) '-' f4gd)) < ERRO;
LOGI1 = LOGI1 'ET' LOGI2;
f5 = 'EXTRAIRE' FLUX4 'H2 ' P1;
ERRO = 'ABS' (1D-12 '*' f5gd) ;
LOGI2 = ('ABS' ((f5 '*' ORIENT) '-' f5gd)) < ERRO;
LOGI1 = LOGI1 'ET' LOGI2;
f6 = 'EXTRAIRE' FLUX4 'O2 ' P1;
ERRO = 'ABS' (1D-12 '*' f6gd) ;
LOGI2 = ('ABS' ((f6 '*' ORIENT) '-' f6gd)) < ERRO;
LOGI1 = LOGI1 'ET' LOGI2;
f7 = 'EXTRAIRE' FLUX4 'H2O ' P1;
ERRO = 'ABS' (1D-12 '*' f7gd) ;
LOGI2 = ('ABS' ((f7 '*' ORIENT) '-' f7gd)) < ERRO;
LOGI1 = LOGI1 'ET' LOGI2;
'SI' ('NON' LOGI1);
'MESSAGE' ;
'MESSAGE' 'OPERATEUR KONV';
'MESSAGE' ('CHAINE' METO);
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'df1 ' ('ABS' (f1gd '-' (ORIENT '*'f1)))
' f1 ' f1);
'MESSAGE' ('CHAINE' 'df2 ' ('ABS' (f2gd '-' (ORIENT '*'f2)))
' f2 ' f2);
'MESSAGE' ('CHAINE' 'df3 ' ('ABS' (f3gd '-' (ORIENT '*'f3)))
' f3 ' f3);
'MESSAGE' ('CHAINE' 'df4 ' ('ABS' (f4gd '-' (ORIENT '*'f4)))
' f4 ' f4);
'MESSAGE' ('CHAINE' 'df5 ' ('ABS' (f5gd '-' (ORIENT '*'f5)))
' f5 ' f5);
'MESSAGE' ('CHAINE' 'df6 ' ('ABS' (f6gd '-' (ORIENT '*'f6)))
' f6 ' f6);
'MESSAGE' ('CHAINE' 'df7 ' ('ABS' (f7gd '-' (ORIENT '*'f7)))
' f7 ' f7);
'ERREUR' 5;
'FINSI' ;
************************************************
***** Splitting des scalaires passifs **********
************************************************
*
LISTINCO = 'MOTS' 'RN' 'RUX' 'RUY' 'RETN' 'H2' 'O2' 'H2O' 'TUTU' 'TATA'
;
TFLUX2 DELTAT = 'KONV' 'VF' 'PERFTEMP' 'FLUX' METO
$DOMTOT PGAZ2 LISTINCO ROF2 VITF2 PF2 YF2 SCAF ;
TFLUX3 = 'EXCO' ( 'MOTS' 'RN' 'RUX' 'RUY' 'RETN' 'H2' 'O2' 'H2O' )
TFLUX2 ( 'MOTS' 'RN' 'RUX' 'RUY' 'RETN' 'H2' 'O2' 'H2O' ) ;
ERRO = 'MAXIMUM' (TFLUX3 '-' TFLUX) 'ABS' ;
'SI' (ERRO > 1.0D-15);
'MESSAGE' 'Probleme dans les transport de scalaires passifs' ;
'ERREUR' 5;
'FINSI' ;
*
**** Test 'KONV'
*
fceltutu = 'EXTRAIRE' TFLUX2 'TUTU' P1;
LOGI1 = ('ABS' ((fceltutu '*' ORIENT) '-' ftutugd)) <
('ABS' (1.0D-12 * fceltutu)) ;
fceltata = 'EXTRAIRE' (TFLUX2) 'TATA' P1;
LOGI2 = ('ABS' ((fceltata '*' ORIENT) '-' ftatagd)) <
('ABS' (1.0D-12 * fceltata)) ;
LOGI1 = LOGI1 'ET' LOGI2;
'SI' ('NON' LOGI1);
'MESSAGE' ;
'MESSAGE' 'OPERATEUR KONV';
'MESSAGE' 'Transport de scalaires passifs' ;
'ERREUR' 5;
'FINSI' ;
****************************************************
****************************************************
******** Fin boucle sur les angles *********
****************************************************
****************************************************
'FIN' BLOC;
'FIN' ;
© Cast3M 2003 - All rights reserved.
Disclaimer