* fichier : couette.dgibi
*----------------------------------------------------------------------
* 1/Ecoulement de Couette engendré par la rotation de 2 cylindres
* concentriques. Comparaison à la solution analytique pour la vitesse.
* 2/Conducion de la chaleur dans un barreau cylindrique (en coordonnées
* cartésiennes. Comparaison à la solution analytique.
* Les calculs sont réalisés successivement avec des éléments
* tri6-(iso-p2) tri7 qua8-(iso-q2) qua9
*----------------------------------------------------------------------
*------------------
* Options générales
*------------------
GRAPH = 'N' ;
COMPLET = FAUX ;
SI ( COMPLET ) ;
nbe=10 ; nbv=20 ;
sinon ;
nbe = 5 ; nbv = 10 ;
finsi ;
ERUT3=0.5 ;
ERTT3=4.E-2;
ERUT7=9.E-2 ;
ERTT7=3.E-2 ;
ERUQ4=7.E-2 ;
ERTQ4=8.E-3 ;
ERUQ9=5.E-2 ;
ERTQ9=2.E-3 ;
DEBPROC TEST ;
ARGU MACRO*MOT GRAPH*MOT EPSS*FLOTTANT
ERRU*FLOTTANT ERRT*FLOTTANT
NBE*ENTIER NBV*ENTIER ;
***********************
*******MAILLAGE********
***********************
**** maillage 1er moitié de cercle ****
* decoupage NBE,NBV
*nbe = 5 ;
*nbv = 10 ;
* Rayon interne
ri = 0.1 ;
drr = 10.;
*rayon externe
re = ri + drr ;
di1=drr/nbe*0.1 ;
di2=drr/nbe*3.;
c1=0. 0. ;
p1o=ri 0. ;
p2o=re 0. ;
tole=2.e-3 ;
*tole=5.e-3 ;
angle=-179.99 ;
angl2= angle/(-2.) ;
p1=p1o tour c1 angl2 ;
p2=p2o tour c1 angl2 ;
q1=p1 tour c1 angle ;
q2=p2 tour c1 angle ;
pp1=p1 c c1 q1 nbv ;
pp2=p2 c c1 q2 nbv ;
mt=pp1 regl dini di1 dfin di2 pp2 ;
entree= cote 4 mt ;
sortie= cote 2 mt ;
axs = entree tour c1 ((-1.)*angl2) ;
elim (axs et mt) tole ;
pp1=inve pp1 ;
sortie=inve sortie ;
cnt1=entree et pp2 et sortie et pp1 ;
couette = mt ;
elim couette tole ;
cntint=pp1 ; elim cntint tole ;
cntext=pp2 ; elim cntext tole ;
*trace couette et cntext et cntint ;
couette=(couette syme DROI c1 (0 100.)) et couette ;
cntint =(cntint syme DROI c1 (0 100.)) et cntint ;
cntext =(cntext syme DROI c1 (0 100.)) et cntext ;
couette = chan couette quaf ;
*trace couette ;
cntint = chan cntint quaf ;
cntext = chan cntext quaf ;
entree = chan entree quaf ;
axs = chan axs quaf ;
elim (couette et cntint et cntext et entree
et axs et sortie ) tole ;
**** definition des contours internes et externes ****
$couette=mode couette 'NAVIER_STOKES' MACRO ;
doma $couette 'IMPR' ;
$cntint = mode cntint 'NAVIER_STOKES' MACRO ;
$cntext = mode cntext 'NAVIER_STOKES' MACRO ;
$entree = mode entree 'NAVIER_STOKES' MACRO ;
$axs = mode axs 'NAVIER_STOKES' MACRO ;
£axs = doma $axs 'TABLE' ;
nbee=(nbno (doma $entree MAILLAGE)) ;
$sortie = mode sortie 'NAVIER_STOKES' MACRO ;
* opti donn 5 ;
************************
*** RESOLUTION ***
************************
* Vitesse tangentielle cylindre interne : u1
* Vitesse tangentielle cylindre externe : u2
u1=-2. ; u2= 2. ;
T1=2. ; T2=1. ;
ris=1./ri ; res=1./re ;
xi= (coor 1 (doma $cntint maillage))*ris ;
yi= (coor 2 (doma $cntint maillage))*ris *(-1.) ;
xe= (coor 1 (doma $cntext maillage))*res ;
ye= (coor 2 (doma $cntext maillage))*res *(-1.) ;
u1x= u1*yi ; u1y= u1*xi ;
u2x= u2*ye ; u2y= u2*xe ;
NU=10. ;
*RV = EQEX $couette OPTI 'EF' 'IMPL' PRESSP0
RV = EQEX $couette OPTI 'EF' 'IMPL'
ZONE $couette OPER LAPN NU INCO 'UN'
ZONE $couette OPER DUDW EPSS INCO 'UN'
'CLIM' 'UN' 'UIMP' cntext u2x 'UN' 'VIMP' cntext u2y
'UN' 'UIMP' cntint u1x 'UN' 'VIMP' cntint u1y ;
mess ' 2 EQEX ' ;
RT = EQEX $couette OPTI 'EF' 'IMPL'
ZONE $couette OPER LAPN NU INCO 'TN'
'CLIM' 'TN' 'TIMP' cntint T1 'TN' 'TIMP' cntext T2 ;
rv.'INCO'=table 'INCO' ;
rv.'INCO'.'UN' = kcht $couette VECT SOMMET (0. 0.) ;
rt.'INCO'=table 'INCO' ;
rt.'INCO'.'TN' = kcht $couette SCAL SOMMET 0. ;
EXEC RV ;
EXEC RT ;
*** Traces ***
si ('EGA' graph 'O' );
un=rv.'INCO'.'UN' ;
tn=rt.'INCO'.'TN' ;
un1=vect un 0.1 ux uy jaune ;
trace un1 couette ;
finsi ;
*** Solutions analytiques ***
****** Comparaison vitesse calculée / vitesse analytique *****
srti=doma $axs 'MAILLAGE' ;
si ( non (ega macro 'MACRO'));
*srtq=chan srti seg2 ;
srtq = srti ;
$axq=mode srtq 'NAVIER_STOKES' MACRO ;
srtv=extr (evol 'CHPO' (ELNO $axq (doma $axq 'VOLUME')) srtq) 'ORDO';
sinon ;
srtv=extr (evol 'CHPO' (ELNO $axs (doma $axs 'VOLUME')) srti) 'ORDO';
finsi ;
uteta=(kops (rv.'INCO'.'UN') 'PSCA' (rv.'INCO'.'UN'))**0.5 ;
evol2 = EVOL 'CHPO' (rv.'INCO'.'UN') UY (srti ) ;
*evol2 = EVOL 'CHPO' uteta 'SCAL' (srti ) ;
r1=evol 'CHPO' (coor 1 srti ) srti ;
r=extr r1 'ORDO' ; $
yu1=extr evol2 'ORDO' ;
evol2= evol 'MANU' 'Rayon' r 'Vitesse' yu1 ;
ric=ri*ri ; rec=re*re ; w1=u1/ri ; w2=u2/re ;
ar=prog nbee* (rec - ric) ;
ue=r*((rec*w2) - (ric*w1)) ;
riec=prog nbee* (ric*rec) ;
ue=ue - ( ( (riec)/r )*(w2-w1) ) ;
ue=ue/ar ;
evol1 = EVOL 'MANU' Rayon r 'Vitesse' ue ;
evolt = evol1 et evol2 ;
TAB1=TABLE;
TAB1.1='MARQ ETOI REGU ' ;
TAB1.2='TIRR REGU ';
TAB1.'TITRE' =TABLE;
TAB1.'TITRE'. 1 =MOT 'Vitesse_théorique ' ;
TAB1.'TITRE'. 2 =MOT 'Vcalculée ';
l1=extr evol1 'ORDO' ; l2=extr evol2 'ORDO' ;
er=(resu ( (l1 - l2 )*(l1 - l2)*srtv ) ) ** 0.5 ;
mess 'Erreur sur la composante azimuthale de la vitesse ' er ;
si ( er > erru ) ; erreur 5 ; finsi ;
si ('EGA' graph 'O' );
DESS evolt 'TITX' 'R (m)' 'TITY' 'V (m/s)' LEGE TAB1 ;
finsi ;
****** Comparaison Teta calculé / Teta analytique *****
kt=(T2-T1)/(Log (Re/Ri)) ;
tl1=prog nbee* T1 ;
rli=prog nbee* Ri ;
teta=tl1 + (kt * (Log(r/Rli))) ;
evl1 = EVOL 'MANU' Rayon r 'Teta_analytique' teta ;
evl2 = EVOL 'CHPO' (rt.'INCO'.'TN') 'SCAL' (srti ) ;
yt1=extr evl2 'ORDO' ;
evl2= evol 'MANU' 'Rayon' r 'Teta_calculé ' yt1;
evlt= evl1 et evl2 ;
TABT=TABLE;
TABT.1='MARQ ETOI REGU ' ;
TABT.2='TIRR REGU ';
TABT.'TITRE'=TABLE;
TABT.'TITRE'. 1=MOT 'Teta_analytique ' ;
TABT.'TITRE'. 2=MOT 'Teta_calculé ';
l1=extr evl1 'ORDO' ; l2=extr evl2 'ORDO' ;
er=(resu ( (l1 - l2 )*(l1 - l2)*srtv ) ) ** 0.5 ;
mess 'Erreur sur la temperature ' er ;
si ( er > errt ) ; erreur 5 ; finsi ;
si ('EGA' graph 'O' );
DESS evlt 'TITX' 'r' 'TITY' 'Teta' LEGE TABT ;
finsi ;
FINPROC RV ;
opti ISOV suli ;
option dime 2 elem tri6 ;
MACRO='MACRO' ;
EPSS=1.E-6 ;
RV=TEST MACRO GRAPH EPSS ERUT3 ERTT3 NBE NBV ;
option dime 2 elem tri6 ;
MACRO='QUAF' ;
EPSS=1.E-12;
RV=TEST MACRO GRAPH EPSS ERUT7 ERTT7 NBE NBV ;
option dime 2 elem qua8 ;
MACRO='MACRO' ;
EPSS=1.E-12 ;
RV=TEST MACRO GRAPH EPSS ERUQ4 ERTQ4 NBV NBV ;
option dime 2 elem qua8 ;
MACRO='QUAF' ;
EPSS=1.E-12 ;
RV=TEST MACRO GRAPH EPSS ERUQ9 ERTQ9 NBV NBV ;
FIN ;
© Cast3M 2003 - All rights reserved.
Disclaimer