Download consistence1_Godunov.dgibi
* fichier : consistence1_Godunov.dgibi
***********************************************************
**** APPROCHE VF "Cell-Centred Formulation" pour la ****
**** solution des ****
**** Equations d'Euler pour un gaz parfait ****
**** ****
**** Gaz monoespece "calorically perfect" ****
**** Consistence, methode GODUNOV ****
**** ****
**** A. BECCANTINI DRN/DMT/SEMT/TTMF NOVEMBRE 1998 ****
***********************************************************
'OPTION' 'DIME' 2 ;
'OPTION' 'ELEM' 'QUA4' ;
'OPTION' 'ISOV' 'SULI' ;
'OPTION' 'ECHO' 0 ;
'OPTION' 'TRAC' 'X';
GRAPH = FAUX ;
* GRAPH = VRAI ;
*
*** Methodes possibles :
*
* 'GODUNOV'
* 'VANLEER'
* 'VLH'
* 'HUSVL'
* 'HUSVLH'
*
METO = 'GODUNOV' ;
************
* MAILLAGE *
************
NX = 5 ;
NY = 1 ;
L = 1.0D0 ;
DX = L/NX/2.0D0 ;
H = DX ;
xD = 0.5D0*L ;
xG = -1.0D0*xD ;
yH = 0.5D0*H ;
yB = -1.0D0*yH ;
A1 = xG yB ;
A2 = 0.0D0 yB ;
A3 = xD yB ;
A4 = xD yH ;
A5 = 0.0D0 yH ;
A6 = xG yH ;
VECTG = XG 0.0D0 ;
VECTD = XD 0.0D0 ;
xBG = xG-DX ;
XBD = xD+DX ;
B1 = xBG yB ;
B2 = xBG yH ;
B3 = xBD yB ;
B4 = xBD yH ;
BAS1 = A1 'DROI' NX A2 ;
BAS2 = A2 'DROI' NX A3 ;
HAU2 = A4 'DROI' NX A5 ;
HAU1 = A5 'DROI' NX A6 ;
LATG = B1 'DROI' NY B2 ;
LAT1 = A1 'DROI' NY A6 ;
LAT12 = A2 'DROI' NY A5 ;
LAT2 = A3 'DROI' NY A4 ;
LATD = B3 'DROI' NY B4 ;
DOM1 = LAT12 'TRANSLATION' NX VECTG ;
DOM2 = LAT12 'TRANSLATION' NX VECTD ;
VECTFG = (-1.0D0*DX) 0.0D0;
VECTFD = DX 0.0D0;
VECTBH = 0.0D0 (yH '-' yB);
VECTHB = 0.0D0 (yB '-' yH);
FRONTG1 = LAT1 'TRANSLATION' 1 VECTFG ;
DCEL1 = FRONTG1 'ET' DOM1;
FRONTG2 = DCEL1 'PLUS' VECTBH;
FRONTG3 = DCEL1 'PLUS' VECTHB;
FRONTG = (FRONTG1 'ET' FRONTG2 'ET' FRONTG3) 'COULEUR' 'ROUG';
'ELIMINATION' 0.0001 FRONTG ;
FRONTD1 = LAT2 'TRANSLATION' 1 VECTFD;
DCEL2 = FRONTD1 'ET' DOM2;
FRONTD2 = DCEL2 'PLUS' VECTBH;
FRONTD3 = DCEL2 'PLUS' VECTHB;
FRONTD = (FRONTD1 'ET' FRONTD2 'ET' FRONTD3) 'COULEUR' 'VERT';
'ELIMINATION' 0.0001 FRONTD ;
*
*** Rotation
*
ANGLE = 30.0D0;
ORIG = 0.0D0 0.0D0;
'MESSAGE' ;
'MESSAGE' ('CHAINE' 'ANGLE = ' ANGLE);
'MESSAGE' ;
DOM1 = DOM1 'TOURNER' ANGLE ORIG;
DOM2 = DOM2 'TOURNER' ANGLE ORIG;
FRONTG = FRONTG 'TOURNER' ANGLE ORIG;
FRONTD = FRONTD 'TOURNER' ANGLE ORIG;
DOMINT = DOM1 'ET' DOM2 ;
'ELIMINATION' DOMINT 1D-6 ;
FRONT = FRONTG ET FRONTD ;
'ELIMINATION' FRONT 1D-6 ;
DOMTOT = DOMINT ET FRONT ;
'ELIMINATION' DOMTOT 1D-6 ;
'SI' GRAPH;
'TRACER' DOMTOT 'TITRE' 'Maillage' ;
'FINSI' ;
$DOMTOT = 'MODE' DOMTOT 'EULER' ;
$DOMINT = 'MODE' DOMINT 'EULER' ;
$DOM1 = 'MODE' DOM1 'EULER' ;
$DOM2 = 'MODE' DOM2 'EULER' ;
$FRONTG = 'MODE' FRONTG 'EULER' ;
$FRONTD = 'MODE' FRONTD 'EULER' ;
$FRONT = 'MODE' FRONT 'EULER' ;
TDOMTOT = 'DOMA' $DOMTOT 'VF' ;
TDOMINT = 'DOMA' $DOMINT 'VF' ;
TDOM1 = 'DOMA' $DOM1 'VF' ;
TDOM2 = 'DOMA' $DOM2 'VF' ;
TFRONTG = 'DOMA' $FRONTG 'VF' ;
TFRONTD = 'DOMA' $FRONTD 'VF' ;
TFRONT = 'DOMA' $FRONT 'VF' ;
MDOMTOT = TDOMTOT . 'QUAF' ;
MDOM1 = TDOM1 . 'QUAF' ;
MDOM2 = TDOM2 . 'QUAF' ;
MDOMINT = TDOMINT . 'QUAF' ;
MFRONTG = TFRONTG . 'QUAF' ;
MFRONTD = TFRONTD . 'QUAF' ;
MFRONT = TFRONT . 'QUAF' ;
'ELIM' (MDOMTOT 'ET' MDOM1 'ET' MDOM2 'ET' MDOMINT 'ET' MFRONTG
'ET' MFRONTD 'ET' MFRONT) 1.E-5 ;
************************************************************
* CONDITIONS INITIALES ET LIMITES. *
************************************************************
gamgd = 1.4D0;
*
*** Etat gauche
*
rog = 1.11;
ung = 21.3;
utg = 11.0;
pg = 1234.7;
rouxg = ((ung '*' ('COS' ANGLE)) '-'
(utg '*' ('SIN' ANGLE))) '*' rog ;
rouyg = ((ung '*' ('SIN' ANGLE)) '+'
(utg '*' ('COS' ANGLE))) '*' rog;
recing = 0.5D0 '*' rog '*' ((ung '*' ung) '+' (utg '*' utg));
retg = (pg '/' (gamgd '-' 1.0)) '+' recing;
*
*** Etat droite
*
rod = 1.11;
und = 21.3;
utd = 11.0;
pd = 1234.7;
rouxd = ((und '*' ('COS' ANGLE)) '-'
(utd '*' ('SIN' ANGLE))) '*' rod;
rouyd = ((und '*' ('SIN' ANGLE)) '+'
(utd '*' ('COS' ANGLE))) '*' rod;
recind = 0.5D0 '*' rod '*' ((und '*' und) '+' (utd '*' utd));
retd = (pd '/' (gamgd '-' 1.0)) '+' recind;
*
*** gamma
*
GAMMA = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' gamgd;
*
*** ro
*
RO_f1 = 'KCHT' $FRONTG 'SCAL' 'CENTRE' rog ;
RO_f2 = 'KCHT' $FRONTD 'SCAL' 'CENTRE' rod ;
RO_f = 'KCHT' $FRONT 'SCAL' 'CENTRE' (RO_f1 'ET' RO_f2) ;
RO_i1 = 'KCHT' $DOM1 'SCAL' 'CENTRE' rog;
RO_i2 = 'KCHT' $DOM2 'SCAL' 'CENTRE' rod;
RO_i = 'KCHT' $DOMINT 'SCAL' 'CENTRE' (RO_i1 'ET' RO_i2);
RN = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' (RO_i 'ET' RO_f) ;
*
*** ro u, ro v
*
GN_f1 = 'KCHT' $FRONTG 'VECT' 'CENTRE' (rouxg rouyg);
GN_f2 = 'KCHT' $FRONTD 'VECT' 'CENTRE' (rouxd rouyd);
GN_f = 'KCHT' $FRONT 'VECT' 'CENTRE' (GN_f1 'ET' GN_f2);
GN_i1 = 'KCHT' $DOM1 'VECT' 'CENTRE' (rouxg rouyg);
GN_i2 = 'KCHT' $DOM2 'VECT' 'CENTRE' (rouxd rouyd);
GN_i = 'KCHT' $DOMINT 'VECT' 'CENTRE' (GN_i1 'ET' GN_i2);
GN = 'KCHT' $DOMTOT 'VECT' 'CENTRE' (GN_i 'ET' GN_f) ;
*
*** ro e
*
RE_f1 = 'KCHT' $FRONTG 'SCAL' 'CENTRE' retg ;
RE_f2 = 'KCHT' $FRONTD 'SCAL' 'CENTRE' retd ;
RE_f = 'KCHT' $FRONT 'SCAL' 'CENTRE' (RE_f1 'ET' RE_f2) ;
RE_i1 = 'KCHT' $DOM1 'SCAL' 'CENTRE' retg ;
RE_i2 = 'KCHT' $DOM2 'SCAL' 'CENTRE' retd ;
RE_i = 'KCHT' $DOMINT 'SCAL' 'CENTRE' (RE_i1 'ET' RE_i2);
RETN = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' (RE_i ET RE_f) ;
********************************************************
*** CREATION DE 'MODE' POUR GRAPHIQUER LE CHAMELEM ***
********************************************************
MOD1 = 'MODELISER' ('DOMA' $DOMTOT 'MAILLAGE') 'THERMIQUE';
*
*** GRAPHIQUE DES C.I.
*
'SI' GRAPH;
*
*** CREATION DE CHAMELEM
*
CHM_RN = 'KCHA' $DOMTOT 'CHAM' RN;
CHM_GN = 'KCHA' $DOMTOT 'CHAM' GN ;
CHM_RETN = 'KCHA' $DOMTOT 'CHAM' RETN;
TRAC CHM_RN MOD1 'TITR' ('CHAINE' 'RHO at t=' 0.0);
TRAC CHM_RETN MOD1 'TITR' ('CHAINE' 'RHOET at t=' 0.0);
TRAC CHM_GN MOD1 'TITR' ('CHAINE' 'RHOU at t=' 0.0);
'FINSI' ;
ZERO =1.0D-12 ;
TFINAL = 100. ;
NITER = 50 ;
SAFFAC = 0.5 ;
*
RN0 = 'COPIER' RN ;
GN0 = 'COPIER' GN ;
RETN0 = 'COPIER' RETN ;
*
TPS = 0. ;
'REPETER' BL1 NITER ;
*
**** Primitive variables
*
VN PN = 'PRIM' 'PERFMONO' RN GN RETN GAMMA ;
ROF VITF PF GAMF = 'PRET' 'PERFMONO' 1 1
$DOMTOT RN VN PN GAMMA;
RESIDU DELTAT = 'KONV' 'VF' 'PERFMONO' 'RESI' METO
$DOMTOT ('MOTS' 'RN' 'RUX' 'RUY' 'RETN')
ROF VITF PF GAMF ;
DT_CON = SAFFAC '*' DELTAT ;
*
**** The time step linked to tps
*
DTTPS = (TFINAL '-' TPS) * (1. '+' ZERO) ;
*
**** Total time step
*
DTMIN = 'MINIMUM' ('PROG' DT_CON DTTPS) ;
*
**** Increment of the variables (convection)
*
RESIDU = DTMIN '*' RESIDU ;
*
DRN = 'REDU' ('EXCO' 'RN' RESIDU 'SCAL') ('DOMA' $DOMINT 'CENTRE');
DGN = 'REDU' ('EXCO' ('MOTS' 'RUX' 'RUY') RESIDU ('MOTS' 'UX' 'UY'))
('DOMA' $DOMINT 'CENTRE') ;
DRETN = 'REDU' ('EXCO' 'RETN' RESIDU 'SCAL')
('DOMA' $DOMINT 'CENTRE') ;
*
TPS = TPS '+' DTMIN ;
RN = RN '+' DRN ;
GN = GN '+' DGN ;
RETN = RETN '+' DRETN ;
*
'SI' (((&BL1 '/' 20) '*' 20) 'EGA' &BL1) ;
'MESSAGE' ('CHAINE' 'ITER =' &BL1 ' TPS =' TPS) ;
'FINSI' ;
'SI' (TPS > TFINAL) ;
'QUITTER' BL1 ;
'FINSI' ;
'FIN' BL1 ;
*
*** GRAPHIQUE DES SOLUTIONS
*
'SI' (GRAPH);
*
*** CREATION DE CHAMELEM
*
CHM_RN = 'KCHA' $DOMTOT 'CHAM' RN ;
CHM_GN = 'KCHA' $DOMTOT 'CHAM' GN ;
CHM_RETN = 'KCHA' $DOMTOT 'CHAM' RETN ;
TRAC CHM_RN MOD1 'TITR' ('CHAINE' 'RHO at t=' TPS);
TRAC CHM_RETN MOD1 'TITR' ('CHAINE' 'RHOET at t=' TPS);
TRAC CHM_GN MOD1 'TITR' ('CHAINE' 'RHOU at t=' TPS);
'FINSI' ;
*
*** Critere de control
*
ERRO = 'MAXIMUM' (RN '-' RN0) 'ABS' ;
'SI' (ERRO '>' 1.0D-6) ;
'ERREUR' 5 ;
'FINSI' ;
'FIN' ;
© Cast3M 2003 - All rights reserved.
Disclaimer