Accueil | Cast3M

* fichier :  flux_wall.dgibi
************************************************************************
************************************************************************
*********************************************************************
**** APPROCHE VF "Cell-Centred Formulation" pour la solution des **** 
**** Equations d'Euler pour un gaz parfait                       ****
**** OPERATEURS PRIM, PRET, KONV                                 ****
**** Gaz monoespece  "calorically perfect"                       ****
**** Flux au mur calcule par Godunov                             ****
****                                                             ****
**** A. BECCANTINI DRN/DMT/SEMT/TTMF    NOVEMBRE  1998           ****
*********************************************************************
 
'OPTION'  'DIME' 2 ;
'OPTION'  'ELEM' QUA4 ;
'OPTION'  'ECHO'  0 ;
'OPTION'  'TRAC' 'PS' ;
 
GRAPH = FAUX ;
* GRAPH = VRAI ;
 
*************************** 
***** DOMAINE SPATIAL  ****
***************************
 
 
A1 = 0.0D0 0.0D0;
A2 = 1.0D0 0.0D0;
A3 = 1.0D0 1.0D0;
A4 = 0.0D0 1.0D0;
 
L12 = A1 'DROIT' 1 A2 ;
L23 = A2 'DROIT' 1 A3 ;
L34 = A3 'DROIT' 1 A4 ;
L41 = A4 'DROIT' 1 A1 ;
 
DOM10 = 'DALL' L12 L23 L34 L41  'PLANE';
 
*
*** Point ou on controlle l'etat 
*
 
P10 = 1.0 0.5 ;
 
*
*** Etat centre
*
 
ro     = 1.0D0 ;
un     = 1.0D0 ;
ut     = 10.0D0 ;
p      = 1.0D0 ;
gam    = 1.4D0 ;
 
roet = (p '/' (gam '-' 1.0)) '+'
      ( 0.5 '*' ro '*' ((un '*' un) '+' (ut '*' ut))) ;
*
*** Etat exact sur le mur
*
 
ps = 2.9266499D0 ;
 
f1ex = 0.0 ;
f2ex = ps ;
f3ex = 0.0 ;
f4ex = 0.0 ;
 
****************************************************
****************************************************
******** Boucle sur les angles             *********
****************************************************
****************************************************
 
DANGLE = 360 '/' 7.15;
ANGLE = 11.3 ;
 
'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;
P1     = P10     'TOURNER' ANGLE ORIG;
 
DOMTOT = DOM1 ;
'ELIM' DOMTOT  1.E-5 ;
 
$DOMTOT = 'MODE' DOMTOT 'EULER'  ;
TDOMTOT = 'DOMA' $DOMTOT 'VF' ;
 
'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  ****
***********************
 
ux = (un '*' ('COS' ANGLE)) '-' (ut '*' ('SIN' ANGLE));
uy = (un '*' ('SIN' ANGLE)) '+' (ut '*' ('COS' ANGLE));
 
GAMMA = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' gam ;
 
RN    = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' ro ;
GN    = 'KCHT' $DOMTOT 'VECT' 'CENTRE' ((ro '*' ux) (ro '*' uy));
EN    = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' roet ;
 
 
***************************
****  L'operateur  PRIM****
***************************
 
VITESSE PRES   = 'PRIM' 'PERFMONO'  RN GN EN GAMMA;
 
***************************
****  L'operateur  PRET****
***************************
 
ORDESP = 1;
ORDTEM = 1;
 
ROF VITF PF GAMF   =  'PRET' 'PERFMONO'  ORDESP ORDTEM
    $DOMTOT RN VITESSE PRES GAMMA;
 
*********************************************************
*** Control des etats sur la surface qui contient P1 ****
*********************************************************
 
GEOP1    = ('DOMA' $DOMTOT 'FACEL') 'ELEM' 'CONTENANT' P1;
GEOP2    = ('DOMA' $DOMTOT 'FACE')  'ELEM' 'CONTENANT' P1;
 
 
ROGEOP1  = 'REDU' ROF GEOP1;
VGEOP1   = 'REDU' VITF GEOP1;
PGEOP1   = 'REDU' PF GEOP1;
GAMGEOP1 = 'REDU' GAMF GEOP1;
REFGEOP1 = 'REDU' VITF GEOP2;
 
********************************
**** 'KONV', TEST GODUNOV   ****
********************************
 
METO = 'GODUNOV' ;
LISTINCO = ('MOTS' 'RN' 'RUNX' 'RUNY' 'RETN') ;
 
CHPFLU DT = 'KONV' 'VF' 'PERFMONO' 'FLUX' METO
   $DOMTOT ROF VITF PF GAMF LISTINCO ;
*
**** Les flux aux interfaces sont dans de
*    CHPOINT FACE
*
 
FLUX1 = 'EXCO' ('EXTRAIRE' LISTINCO 1) CHPFLU ;
FLUX2X = 'EXCO' ('EXTRAIRE' LISTINCO 2) CHPFLU ;
FLUX2Y = 'EXCO' ('EXTRAIRE' LISTINCO 3) CHPFLU ;
FLUX3 = 'EXCO' ('EXTRAIRE' LISTINCO 4) CHPFLU ;
 
FLUX2N = (FLUX2X '*' ('COS' ANGLE)) '+' (FLUX2Y * ('SIN' ANGLE));
FLUX2T = (FLUX2Y '*' ('COS' ANGLE)) '-' (FLUX2X * ('SIN' ANGLE));
 
*
**** Test 'KONV'
*
 
ERRO = 1D-6 ;
 
f1 = 'EXTRAIRE' FLUX1 'SCAL' P1;
ERROF1 = 'ABS' (f1 '-' f1ex);
LOGI1 = ERROF1 > ERRO;
 
'SI' LOGI1 ;
  'MESSAGE' ;
  'MESSAGE' ('CHAINE' 'ERR F1 ' ERROF1);
  'MESSAGE' ;
  'ERREUR' 5 ;
'FINSI' ;  
 
f2 = 'EXTRAIRE' FLUX2N 'SCAL' P1;
ERROF2 = 'ABS' (f2 '-' f2ex);
LOGI1 = ERROF2 > ERRO;
 
'SI' LOGI1 ;
  'MESSAGE' ;
  'MESSAGE' ('CHAINE' 'ERR F2 ' ERROF2);
  'MESSAGE' ;
  'ERREUR' 5 ;
'FINSI' ;  
 
f3 = 'EXTRAIRE' FLUX2T 'SCAL' P1;
ERROF3 = 'ABS' (f3 '-' f3ex); 
LOGI1 =  ERROF3 > ERRO;
 
'SI' LOGI1 ;
  'MESSAGE' ;
  'MESSAGE' ('CHAINE' 'ERR F3 ' ERROF3);
  'MESSAGE' ;
  'ERREUR' 5 ;
'FINSI' ;  
 
f4 = 'EXTRAIRE' FLUX3 'SCAL' P1;
ERROF4 = 'ABS' (f4 '-' f4ex);
LOGI1 = ERROF4 > ERRO;
 
'SI' LOGI1 ;
  'MESSAGE' ;
  'MESSAGE' ('CHAINE' 'ERR F4 ' ERROF4);
  'MESSAGE' ;
  'ERREUR' 5 ;
'FINSI' ;  
 
 
'FIN' BLOC ;
 
'FIN' ;