* fichier : fcourant.dgibi
'OPTI' 'ECHO' 0 ;
************************************************************************
* NOM : FCOURANT
* DESCRIPTION : Teste la procédure calculant la fonction de courant
* sur un poiseuille en 2D et 2D axi
*
*
*
* LANGAGE : GIBIANE-CAST3M
* AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
* mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION : v1, 21/12/2007, version initiale
* HISTORIQUE : v1, 21/12/2007, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* Prière de PRENDRE LE TEMPS de compléter les commentaires
* en cas de modification de ce sous-programme afin de faciliter
* la maintenance !
************************************************************************
'SAUTER' 2 'LIGNE' ;
'MESSAGE' ' Execution de fcourant.dgibi' ;
'SAUTER' 2 'LIGNE' ;
*
interact = FAUX ;
graph = FAUX ;
*
'OPTION' 'DIME' 2 'ELEM' 'TRI6' ;
*
'SI' ('NON' interact) ;
'OPTION' 'TRAC' 'PS' ;
'SINON' ;
'OPTION' 'TRAC' 'X' ;
'FINSI' ;
*
lok = VRAI ;
*
'REPETER' iaxi 2 ;
iiaxi = &iaxi ;
'SAUTER' 1 'LIGNE' ;
'SI' ('EGA' iiaxi 1) ;
'OPTI' 'MODE' 'PLAN' ;
'MESSAGE' 'Mode plan ' ;
'SINON' ;
'OPTI' 'MODE' 'AXIS' ;
'MESSAGE' 'Mode axi ' ;
'FINSI' ;
'SAUTER' 1 'LIGNE' ;
*
laxi = ('EGA' ('VALEUR' 'MODE') 'AXIS') ;
*
* Maillage
*
a = 1.2 ; b = 2.5 ;
nx = 5 ; ny = 7 ;
*nx = 10 ; ny = 14 ;
dpdz = -1.5 ; mu = 0.9 ;
vtoim = 0. ('*' dpdz b) ;
disv = 'QUAF' ;
disp = 'CENTREP1' ;
*
pA = 0. 0. ; pB = a 0. ; pC = a b ; pD = 0. b ;
bas = 'DROIT' nx pA pB ; dro = 'DROIT' ny pB pC ;
hau = 'DROIT' nx pC pD ; gau = 'DROIT' ny pD pA ;
*
cmt = bas 'ET' dro 'ET' hau 'ET' gau ;
mt = 'SURFACE' cmt ;
_mt = 'CHANGER' mt 'QUAF' ;
*
$mt = 'MODE' _mt 'NAVIER_STOKES' disv ;
$hau = 'MODE' hau 'NAVIER_STOKES' disv ;
*
mt = 'DOMA' $mt 'MAILLAGE' ;
*
* Solution exacte 2D
* u_x = 0
* u_y = vmax (1 - x^2/a^2)
* vmax = -dp/dz (a^2 / 2 mu)
* p = dp/dz y
* psi = vmax x (1 - x^2/3a^2)
*
'SI' ('NON' laxi) ;
vmax = '*' ('*' dpdz -1.) ('/' ('*' a a) ('*' mu 2)) ;
x = 'COORDONNEE' 1 mt ;
y = 'COORDONNEE' 2 ('DOMA' $mt disp) ;
cy = '+' ('*' ('/' ('*' x x) ('*' a a)) -1.) 1. ;
unx = 'NOMC' 'UY' ('*' vmax cy) ;
pnx = '*' dpdz y ;
cpsi = '+' ('*' ('/' ('*' x x) ('*' ('*' a a) 3)) -1.) 1. ;
psix = '*' vmax ('*' x cpsi) ;
'SINON' ;
*
* Solution exacte 2D axi
* u_r = 0
* u_z = vmax (1 - r^2/a^2)
* vmax = -dp/dz (a^2 / 4 mu)
* p = dp/dz y
* psi = vmax x (1 - x^2/3a^2)
*
vmax = '*' ('*' dpdz -1.) ('/' ('*' a a) ('*' mu 4)) ;
r = 'COORDONNEE' 1 mt ;
z = 'COORDONNEE' 2 ('DOMA' $mt disp) ;
cz = '+' ('*' ('/' ('*' r r) ('*' a a)) -1.) 1. ;
unx = 'NOMC' 'UY' ('*' vmax cz) ;
pnx = '*' dpdz z ;
cpsi = '+' ('*' ('/' ('*' r r) ('*' ('*' a a) 2)) -1.) 1. ;
psix = '*' ('*' vmax PI) ('*' ('*' r r) cpsi) ;
'FINSI' ;
*
* Poiseuille
*
rv = 'EQEX' 'ITMA' 1 'ALFA' 1. 'NITER' 1
'OPTI' 'EF' 'IMPL' 'CENTREE' disp
'ZONE' $mt 'OPER' 'NS' 1. 'UN' mu 'INCO' 'UN'
'ZONE' $mt 'OPER' 'KBBT' 1. 'INCO' 'UN' 'PN'
'ZONE' $hau 'OPER' 'TOIM' vtoim 'INCO' 'UN'
'CLIM' 'UN' 'UIMP' dro 0. 'UN' 'VIMP' dro 0.
'CLIM' 'UN' 'UIMP' gau 0.
;
rv . 'INCO' = 'TABLE' 'INCO' ;
RV . 'INCO' . 'UN' = 'KCHT' $mt 'VECT' 'SOMMET' (0. 0.) ;
RV . 'INCO' . 'PN' = 'KCHT' $mt 'SCAL' disp 0. ;
*
EXEC rv ;
*
un = rv . 'INCO' . 'UN' ;
pn = rv . 'INCO' . 'PN' ;
clpsi = 'MANUEL' 'CHPO' pA 1 'SCAL' 0. ;
psi = FCOURANT _mt disv un clpsi ;
*
'SI' graph ;
titg = 'CHAINE' ' mod = ' ('VALEUR' 'MODE') ;
vun = 'VECT' un 'UX' 'UY' 'JAUN' ;
tit = 'CHAINE' 'Vitesse' titg ;
'TRACER' vun mt cmt 'TITR' tit ;
ppn = 'ELNO' $mt pn disp ;
tit = 'CHAINE' 'Pression' titg ;
'TRACER' ppn mt 'TITR' tit ;
tit = 'CHAINE' 'Fonction de courant' titg ;
'TRACER' psi mt 'TITR' tit ;
'FINSI' ;
*
erru = 'MAXIMUM' ('-' un unx) 'ABS' ;
errp = 'MAXIMUM' ('-' pn pnx) 'ABS' ;
errf = 'MAXIMUM' ('-' psi psix) 'ABS' ;
*
'MESSAGE' ('CHAINE' 'erru = ' erru) ;
'MESSAGE' ('CHAINE' 'errp = ' errp) ;
'MESSAGE' ('CHAINE' 'errf = ' errf) ;
*
* le 21/12/2007 sur ma machine :
* Mode plan
* erru = 4.21885E-15
* errp = 3.15303E-14
* errf = 6.36767E-04
* Mode axi
* erru = 2.66454E-15
* errp = 2.62013E-14
* errf = 5.45703E-03
*
tolup = 1.D-12 ;
'SI' laxi ;
tolf = 6.D-3 ;
'SINON' ;
tolf = 7.D-4 ;
'FINSI' ;
*
tst = ('<' erru tolup) ;
'SI' ('NON' tst) ;
chmes = 'CHAINE' '!!! Erreur, on na pas erru < ' tolup ;
'MESSAGE' chmes ;
'FINSI' ;
lok = lok 'ET' tst ;
*
tst = ('<' errp tolup) ;
'SI' ('NON' tst) ;
chmes = 'CHAINE' '!!! Erreur, on na pas errp < ' tolup ;
'MESSAGE' chmes ;
'FINSI' ;
lok = lok 'ET' tst ;
*
tst = ('<' errf tolf) ;
'SI' ('NON' tst) ;
chmes = 'CHAINE' '!!! Erreur, on na pas errf < ' tolf ;
'MESSAGE' chmes ;
'FINSI' ;
lok = lok 'ET' tst ;
*
'FIN' iaxi ;
*
'SAUTER' 2 'LIGNE' ;
'SI' lok ;
'MESSAGE' 'Tout sest bien passe' ;
'SINON' ;
'MESSAGE' 'Il y a eu des erreurs' ;
'FINSI' ;
'SAUTER' 2 'LIGNE' ;
*
'SI' interact ;
'OPTION' 'DONN' 5 'ECHO' 1 ;
'FINSI' ;
'SI' ('NON' lok) ;
'ERREUR' 5 ;
'FINSI' ;
*
* End of dgibi file FCOURANT
*
'FIN' ;
© Cast3M 2003 - All rights reserved.
Disclaimer