Annexe A : codes de calcul

************************
* Tracé en Postscript Couleur *
************************

OPTION 'TRAC' PSC;

*****************
* Choix de la durée *
*****************

COMPLET = FAUX;
SI ( COMPLET ) ;
d1= 0.5 ; d2= 0.02 ;
nbit = 25000 ;
err1 = 7.5E-3 ;
SINON ;
d1 = 1.0 ; d2 = 0.1 ;
nbit = 5000 ;
err1 = 7.E-2 ;
FINSI ;
GRAPH = 'N' ;


******************************
* Calcul de variations log[Un-(Un-1)] *
******************************

DEBPROC CALCUL ;
ARGU RVX*'TABLE' ;

EXPORT RVX;

RV = RVX.'EQEX' ;

DD = RV.PASDETPS.'NUPASDT' ;
NN = DD/50 ;
LO = (DD-(50*NN)) EGA 0 ;
SI ( LO ) ;

UN = RV.INCO.'UN' ;
UNM1 = RV.INCO.'UNM1' ;

vch= kcht $MT scal sommet (exco 'UX' RV.INCO.'UN') ;

unx= kcht $MT scal sommet (exco 'UX' un) ;
unm1x = kcht $MT scal sommet (exco 'UX' unm1) ;
uny= kcht $MT scal sommet (exco 'UY' un) ;
unm1y = kcht $MT scal sommet (exco 'UY' unm1) ;

ERRX = KOPS unx '-' unm1x ;
ERRY = KOPS uny '-' unm1y ;
ELIX = MAXI ERRX 'ABS' ;
ELIY = MAXI ERRY 'ABS' ;
ELIX = (LOG (ELIX + 1.0E-10))/(LOG 10.) ;
ELIY = (LOG (ELIY + 1.0E-10))/(LOG 10.) ;
MESSAGE 'ITER ' RV.PASDETPS.'NUPASDT' ' ERREUR LINF ' ELIX ELIY ;
RV.INCO.'UNM1' = KCHT $MT 'VECT' 'SOMMET' (RV.INCO.'UN') ;
IT = PROG RV.PASDETPS.'NUPASDT' ;
ER = PROG ELIY ;
RV.INCO.'IT' = (RV.INCO.'IT') ET IT ;
RV.INCO.'ER' = (RV.INCO.'ER') ET ER ;
FINSI ;

FINPROC;

**********************
* exportation des données *
**********************

DEBPROC EXPORT ;
ARGU RVX*'TABLE' ;

RV=RVX.'EQEX';

DD = RV.PASDETPS.'NUPASDT' ;
NN = DD/100 ;
LO = (DD-(100*NN)) EGA 0 ;
SI ( LO ) ;

UN = RV.INCO.'UN' ;

vch= kcht $MT scal sommet (exco 'UX' RV.INCO.'UN');

*Température*
t=RV.INCO.'TN';
nomfic = chaine 'Temp' DD '.ps';
OPTION 'FTRAC' nomfic;
trac t mt;
MESSAGE 'OK pour ' nomfic;

*Vitesses*
unch = vect UN 0.3 UX UY ROUGE ;
nomfic = chaine 'Vit' DD '.ps';
OPTION 'FTRAC' nomfic;
TRACE UNCH MT ;
MESSAGE 'OK pour ' nomfic;

*Coupe vitesses Axiale*
nomfic = chaine 'VitAxiale' DD '.ps';
OPTION 'FTRAC' nomfic;
l00 = droit p0 80 p00;
UU = exco UN UX;
UU = int_comp mt UU l00;
dessin (evol CHPO UU l00) 'YBOR' 0. 5.;
MESSAGE 'OK pour ' nomfic;

FINSI ;

FINPROC;

*************
* MAILLAGE *
*************

option dime 2 elem qua8 ;

p1= 0 -5;
p2= 20 -15;
p3= 20 -3;
p4 =20 3;
p5=20 15;
p6=0 5;
p7=0 1;
p8=0 -1;

p0 = 0 0;
p00=20 0;

ns=5;
nd=50;
d1=0.2;
d2=2*d1;
d3=3*d1;
d4=2*d3;

l12 = droi p1 nd p2 dini d2 dfin d4;
l23 = droi p2 p3 dini d4 dfin d3;
l34 = droi p3 p4 dini d3 dfin d3;
l45 = droi p4 p5 dini d3 dfin d4;
l56 = droi p5 nd p6 dini d4 dfin d2;
l67 = droi p6 p7 dini d2 dfin d1;
l78 = droi p7 p8 dini d1 dfin d1;
l81 = droi p8 p1 dini d1 dfin d2;


entree=l78;
sortie=(l23 et l34 et l45);
mt= l12 (l23 et l34 et l45) l56 (l67 et l78 et l81) 'DALL';
tass mt;
orienter mt;


$mt ='DOMA' mt 1.E-7 'MACRO';
$INLET = 'DOMA' entree 'MACRO' 'INCL' $mt 1.E-7;
$OUTLET = 'DOMA' sortie 'MACRO' 'INCL' $mt 1.E-7;
$l12 = 'DOMA' l12 'MACRO' 'INCL' $mt 1.E-7;
$l56 = 'DOMA' l56 'MACRO' 'INCL' $mt 1.E-7;

**************
* Profil d'entree *
**************

Umax = 4 ;
Yin = 'COOR' 2 $INLET.MAILLAGE;;
Ymax = 'MAXI' Yin;
Ymin = 'MINI' Yin;
Uin = KOPS (KOPS Yin '-' Ymax) '*' (KOPS Yin '-' Ymin);
Uin = KOPS (Uin) '*' (-1.*Umax/((Ymax-Ymin)*(Ymax-Ymin)));
Uin = NOMC 'UX' Uin 'NATU' 'DISCRET';
Vin = kcht $INLET 'SCAL' 'SOMMET' 0.;
Vin = NOMC 'UY' Vin 'NATU' 'DISCRET';

********
* Calcul *
********
Nu = 1./50. ;
tos = 0. 0. ;
gb= -0.1 0.;
ALF=1./50.;
TREF=0.;

RV = EQEX $mt 'ITMA' nbit 'ALFA' 0.8
'OPTI' 'SUPG'
'ZONE' $mt 'OPER' 'NS' Nu gb 'TN' 0. 'INCO' 'UN'
'ZONE' $mt 'OPER' 'TSCAL' ALF 'UN' 0. 'INCO' 'TN'
'ZONE' $OUTLET 'OPER' 'TOIM' tos 'INCO' 'UN';

RV = EQEX RV 'ZONE' $MT 'OPER' CALCUL;

*Conditions sur le sol*

RV = EQEX RV 'CLIM' 'UN' 'UIMP' $INLET.MAILLAGE Uin ;
RV = EQEX RV 'CLIM' 'UN' 'VIMP' $INLET.MAILLAGE Vin ;
RV = EQEX RV 'CLIM' 'TN' 'TIMP' $INLET.MAILLAGE 1.;

RV = EQEX RV 'CLIM' 'UN' 'UIMP' l67 0.;
RV = EQEX RV 'CLIM' 'UN' 'VIMP' l67 0.;

RV = EQEX RV 'CLIM' 'UN' 'UIMP' l81 0.;
RV = EQEX RV 'CLIM' 'UN' 'VIMP' l81 0.;

RV = EQEX RV 'CLIM' 'TN' 'TIMP' l67 0.;
RV = EQEX RV 'CLIM' 'TN' 'TIMP' l81 0.;

*Conditions sur les cotés, soit TOIM, soit UIMP,VIMP*

RV = EQEX RV 'ZONE' $l12 'OPER' 'TOIM' tos 'INCO' 'UN';
RV = EQEX RV 'ZONE' $l56 'OPER' 'TOIM' tos 'INCO' 'UN';

*RV = EQEX RV 'CLIM' 'UN' 'UIMP' l12 0.;
*RV = EQEX RV 'CLIM' 'UN' 'VIMP' l12 0.;

*RV = EQEX RV 'CLIM' 'UN' 'UIMP' l56 0.;
*RV = EQEX RV 'CLIM' 'UN' 'VIMP' l56 0.;

*******************************

RVP = EQPR $mt KTYPI 1
ZONE $mt 'OPER' 'PRESSION' 0.;

RV.INCO = TABLE INCO ;
RV.INCO.'UN' = KCHT $mt VECT SOMMET (1.E-3 1.E-3 ) ;
RV.INCO.'TN' = KCHT $mt SCAL SOMMET (1.E-3 ) ;
RV.INCO.'UNM1' = KCHT $mt 'VECT' 'SOMMET' (1.E-5 1.E-5 ) ;
RV.INCO.'IT' = PROG 1 ;
RV.INCO.'ER' = PROG 0. ;
RV.PRESSION = RVP ;

*historique sur une ligne*

lh = mt POIN 'DROIT' p0 p00 1.e-3;
his = khis 'UN' 1 lh;
rv.'HIST' = his;

*******************************
EXEC RV ;
*******************************

**********
* Resultats *
**********

* Vitesse sur l'axe
*l00 = droit p0 50 p00;
*UU = exco rv.inco.'UN' UX;
*UU = int_comp mt UU l00;
*dessin (evol CHPO UU l00);


* Champ de vitesses
*unch = vect RV.INCO.'UN' 0.4 UX UY JAUNE ;
*TRACE UNCH MT;

* Variations de la vitesse sur un domaine
*dessin his.'TABD' his.'1UN';

* Température
*t = kcht $MT scal sommet rv.inco.'TN';
*trac t mt;


Dans ce cas, on utilise une version légèrement différente pour le maillage dont on ne prend qu’un coté par rapport à l’axe central.

************************
* Tracé en Postscript Couleur *
************************

OPTION 'TRAC' PSC;

M = 0;
Nimages = 50;

****************
* Choix de la durée *
****************

COMPLET = FAUX;
SI ( COMPLET ) ;
d1= 0.5 ; d2= 0.02 ;
nbit = 25000 ;
err1 = 7.5E-3 ;
SINON ;
d1 = 1.0 ; d2 = 0.1 ;
nbit = 5000 ;
err1 = 7.E-2 ;
FINSI ;
GRAPH = 'N' ;


******************************
* Calcul de variations log[Un-(Un-1)] *
******************************

DEBPROC CALCUL ;
ARGU RVX*'TABLE' ;

RV = RVX.'EQEX' ;

DD = RV.PASDETPS.'NUPASDT' ;
NN = DD/50 ;
LO = (DD-(50*NN)) EGA 0 ;
SI ( LO ) ;

UN = RV.INCO.'UN' ;
UNM1 = RV.INCO.'UNM1' ;

vch= kcht $MT scal sommet (exco 'UX' RV.INCO.'UN') ;

unx= kcht $MT scal sommet (exco 'UX' un) ;
unm1x = kcht $MT scal sommet (exco 'UX' unm1) ;
uny= kcht $MT scal sommet (exco 'UY' un) ;
unm1y = kcht $MT scal sommet (exco 'UY' unm1) ;

ERRX = KOPS unx '-' unm1x ;
ERRY = KOPS uny '-' unm1y ;
ELIX = MAXI ERRX 'ABS' ;
ELIY = MAXI ERRY 'ABS' ;
ELIX = (LOG (ELIX + 1.0E-10))/(LOG 10.) ;
ELIY = (LOG (ELIY + 1.0E-10))/(LOG 10.) ;
MESSAGE 'ITER ' RV.PASDETPS.'NUPASDT' ' ERREUR LINF ' ELIX ELIY ;
RV.INCO.'UNM1' = KCHT $MT 'VECT' 'SOMMET' (RV.INCO.'UN') ;
IT = PROG RV.PASDETPS.'NUPASDT' ;
ER = PROG ELIY ;
RV.INCO.'IT' = (RV.INCO.'IT') ET IT ;
RV.INCO.'ER' = (RV.INCO.'ER') ET ER ;
FINSI ;

FINPROC;


*************
* MAILLAGE *
*************

option dime 2 elem qua8 ;

p1= 0 0;
p2= 40 0;
p3= 40 3;
p4 =40 15;
p5=0 5;
p6=0 1;

nd=50;
d1=0.2;
d2=2*d1;
d3=3*d1;
d4=2*d3;

l12 = droi p1 nd p2 dini d1 dfin d3 ;
l23 = droi p2 p3 dini d3 dfin d3;
l34 = droi p3 p4 dini d3 dfin d4;
l45 = droi p4 nd p5 dini d4 dfin d2;
l56 = droi p5 p6 dini d2 dfin d1;
l61 = droi p6 p1 dini d1 dfin d1;


entree = l61;
sortie = (l23 et l34);
mt= l12 (l23 et l34) l45 (l56 et l61) 'DALL';

tass mt;
orienter mt;


$mt ='DOMA' mt 1.E-7 'MACRO';
$INLET = 'DOMA' entree 'MACRO' 'INCL' $mt 1.E-7;
$OUTLET = 'DOMA' sortie 'MACRO' 'INCL' $mt 1.E-7;
$cote = 'DOMA' l45 'MACRO' 'INCL' $mt 1.E-7;

**************
* Profil d'entree *
**************

Umax = 4 ;
Yin = 'COOR' 2 $INLET.MAILLAGE;;
Ymax = 'MAXI' Yin;
Ymin = -1*Ymax;
Uin = KOPS (KOPS Yin '-' Ymax) '*' (KOPS Yin '-' Ymin);
Uin = KOPS (Uin) '*' (-1.*Umax/((Ymax-Ymin)*(Ymax-Ymin)));
Uin = NOMC 'UX' Uin 'NATU' 'DISCRET';
Vin = kcht $INLET 'SCAL' 'SOMMET' 0.;
Vin = NOMC 'UY' Vin 'NATU' 'DISCRET';

Psi1n = KOPS (KOPS (KOPS Yin '*' Yin) '*' Yin) '*' (-1./3.);
Psi2n = KOPS (KOPS Yin '*' Yin) '*' ((Ymax+Ymin)/2.);
Psi3n = KOPS (Yin) '*' (-1*Ymax*Ymin);
Psin = KOPS (Psi1n) '+' (Psi2n);
Psin = KOPS (Psin) '+' (Psi3n);
Psin = KOPS (Psin) '*' (Umax/((Ymax-Ymin)*(Ymax-Ymin)));
Psin = NOMC 'PSI' Psin 'NATU' 'DISCRET';

********
* Calcul *
********
Nu = 1./50. ;
tos = 0. 0. ;
gb= -1. 0.;
ALF=1./50.;
TREF=0.;

RV = EQEX $mt 'ITMA' nbit 'ALFA' 0.8
'OPTI' 'SUPG'
'ZONE' $mt 'OPER' 'NS' Nu gb 'TN' 0. 'INCO' 'UN'
'ZONE' $mt 'OPER' 'TSCAL' ALF 'UN' 0. 'INCO' 'TN'
'ZONE' $OUTLET 'OPER' 'TOIM' tos 'INCO' 'UN';

RV = EQEX RV 'ZONE' $MT 'OPER' CALCUL;

*Conditions sur le sol*

RV = EQEX RV 'CLIM' 'UN' 'UIMP' $INLET.MAILLAGE Uin ;
RV = EQEX RV 'CLIM' 'UN' 'VIMP' $INLET.MAILLAGE Vin ;
RV = EQEX RV 'CLIM' 'TN' 'TIMP' $INLET.MAILLAGE 1.;

RV = EQEX RV 'CLIM' 'UN' 'UIMP' l56 0.;
RV = EQEX RV 'CLIM' 'UN' 'VIMP' l56 0.;

RV = EQEX RV 'CLIM' 'TN' 'TIMP' l56 0.;

*Conditions sur les cotés*

*RV = EQEX RV 'ZONE' $cote 'OPER' 'TOIM' tos 'INCO' 'UN';
RV = EQEX RV 'CLIM' 'UN' 'UIMP' l45 0.;
RV = EQEX RV 'CLIM' 'UN' 'VIMP' l45 0.;

*Symétrisation sur l12 pour éviter Coanda*

RV = EQEX RV 'CLIM' 'UN' 'VIMP' l12 0.;

*******************************

RVP = EQPR $mt KTYPI 1
ZONE $mt 'OPER' 'PRESSION' 0.;

RV.INCO = TABLE INCO ;
RV.INCO.'UN' = KCHT $mt VECT SOMMET (1.E-3 1.E-3 ) ;
RV.INCO.'TN' = KCHT $mt SCAL SOMMET (1.E-3 ) ;
RV.INCO.'UNM1' = KCHT $mt 'VECT' 'SOMMET' (1.E-5 1.E-5 ) ;
RV.INCO.'IT' = PROG 1 ;
RV.INCO.'ER' = PROG 0. ;
RV.PRESSION = RVP ;



REPETER BLOC Nimages;

*******************************
EXEC RV ;
*******************************

**********
* Resultats *
**********
M = M + nbit ;

* Champ de vitesses
unch = vect RV.INCO.'UN' 0.4 UX UY ROUGE ;
nomfic = chaine 'Vit' M '.ps';
OPTION 'FTRAC' nomfic;
TRACE UNCH MT;
MESSAGE 'OK pour ' nomfic;

****************************************************************************
* ON CALCULE LA FONCTION DE COURANT NUMERIQUEMENT A PARTIR DE LA *
* DONNEE DU CHAMP DE VITESSE lignedecourant.dgbi *
****************************************************************************

CNT = CONT mt ;
UUU=rv.inco.'UN';
SW = (-1.) * (KOPS UUU 'ROT' $mt) ;
RK = EQEX $mt 'OPTI' 'EF' 'IMPL'
ZONE $mt OPER LAPN 1.0 INCO 'PSI'
ZONE $mt OPER FIMP sw INCO 'PSI';
RK = EQEX RK 'CLIM' 'PSI' 'TIMP' (l45 et l56) (2/3.);
RK = EQEX RK 'CLIM' 'PSI' 'TIMP' $INLET.MAILLAGE Psin;
RK = EQEX RK 'CLIM' 'PSI' 'TIMP' l12 (0.);

RK.'INCO' = TABLE INCO ;
RK.'INCO'.'PSI' = KCHT $mt SCAL SOMMET (2./3.) ;
EXEC RK ;

PSI = RK.'INCO'.'PSI' ;

nomfic = chaine 'Psi' M '.ps';
OPTION 'FTRAC' nomfic;
TRACE PSI mt CNT 28;
MESSAGE 'OK pour ' nomfic;
FIN BLOC;

*Coupe vitesses Axiale
nomfic = chaine 'VitAxiale' '.ps';
OPTION 'FTRAC' nomfic;
l00 = droit p1 80 p2;
UU = exco RV.INCO.UN UX;
UU = int_comp mt UU l00;
dessin (evol CHPO UU l00) 'YBOR' 0. 6.;

i= 0;
M = 20;
UU = exco RV.INCO.UN UX;
REPETER BLOU M;
nomfic = chaine 'ProfilVit' i '.ps';
OPTION 'FTRAC' nomfic;
pA = 2*i 0;
pB = 2*i 15;
l000 = droit pA 80 pB;
UUU = int_comp mt UU l000;
dessin (evol CHPO UUU l000) 'XBOR' 0 6 'YBOR' 0. 6.;
MESSAGE 'OK pour ' nomfic;
i = i+1;
FIN BLOU;


UU = exco RV.INCO.UN UX;
UU = int_comp mt UU l00;
list UU;

**** fin ****