* fichier : Hpressu.dgibi
*---------------------------------------------------------
* Pressurisation d'une enceinte type Phébus
*
* Le maillage correspond à une enceinte cylindrique
* d'environ 10 m3 avec un mur en contact avec la
* paroi verticale (10 cm d'acier)
* Tout le volume est initialement a 1bar et 40oC
* et la température du mur est mise à 60oC
* au début du calcul. On calcule la pressurisation
* de cette enceinte sur 50 secondes en injectant un débit
* de 50g/s de vapeur à 150oC. Ce test (un peu long)
* verifie le demarrage de la condensation a 20"
* verifie l'evolution moyenne de la temprature gaz
* verifie la pression max a 50"
* verifie la vitesse max a 50" (Convection naturelle)
* verifie la masse d'eau a 50"
* Auteurs : E. Studer, J.P. Magnaud Novembre 1999
* revisite Mars 2001
*--------------------------------------------------------
COMPLET = FAUX ;
*COMPLET = VRAI ;
GRAPH = VRAI ;
GRAPH = FAUX ;
'SI' COMPLET ;
nbit=100;
DT0 = 1. ;
n1 = 1 ; n2 = 4 ; n3 = 4 ;
n4 = 8 ; nn = 2 ;
'SINON' ;
nbit= 5 ;
DT0 = 10. ;
n1 = 1 ; n2 = 2 ; n3 = 4 ;
n4 = 4 ; nn = 1 ;
'FINSI' ;
*--------------------------------------------------------
* Definition du maillage de l'enceinte cylindrique
*
'OPTI' 'DIME' 3 'ELEM' 'CU20' ;
'OPTI' 'TRAC' 'X' ;
ri = 1.052 ; sp = 0.10 ; hc = 4.163 ;
epsi = 1.000e-2 ; ;
p0 = 0.000 0.000 0.000 ;
px = -1000.000 0.000 0.000 ;
py = 0.000 -1000.000 0.000 ;
pz = 0.000 0.000 1000.000 ;
cd = 0.000 0.000 -20.000 ;
ph0 = 0.000 0.000 hc ;
phx = ri 0.000 hc ;
phy = 0.000 ri hc ;
fg1 = 0.25 ;
fg2 = fg1 * (2.0 ** 0.5) / 2. ;
p1 = (ri*fg1) 0.000 0.000 ;
p2 = (ri*fg2) (ri*fg2) 0.000 ;
p3 = 0.000 (ri*fg1) 0.000 ;
p4 = ri 0.000 0.000 ;
p5 = 0.000 ri 0.000 ;
p6 = (ri+sp) 0.000 0.000 ;
p7 = 0.000 (ri+sp) 0.000 ;
* Hauteur de l'enceinte
h1 = 4.163 ;
* Vecteur de translation
v1 = 0. 0. h1 ;
l1 = 'DROI' p0 p1 n1 ;
l2 = 'DROI' p1 p2 n1 ;
l3 = 'DROI' p2 p3 n1 ;
l4 = 'DROI' p3 p0 n1 ;
l5 = 'CERC' p4 p0 p5 (2*n1) ;
l6 = 'CERC' p6 p0 p7 (2*n1) ;
basf0= 'DALL' l1 l2 l3 l4 'PLAN' ;
basf1=('REGL' (l2 'ET' l3) l5 n2) ;
l44= cote 2 basf1;
ax4= (inve l4) et l44 ;
l11= cote 4 basf1;
ax1= l11 et (inve l1) ;
basf = basf0 'ET' ('REGL' (l2 'ET' l3) l5 n2) ;
'ELIM' basf epsi ;
basf = basf 'ET' ('SYME' basf 'DROI' p0 p3) ;
ax11 = ('SYME' ax1 'DROI' p0 p3) 'ET' (inve ax1) ;
'ELIM' basf epsi ;
basf = basf 'ET' ('SYME' basf 'DROI' p0 p1) ;
ax44 = (inve ax4) 'ET' ('SYME' ax4 'DROI' p0 p4) ;
'ELIM' basf epsi ;
basm = 'REGL' l5 l6 n3 ;
basm = basm 'ET' ('SYME' basm 'DROI' p0 p3) ;
'ELIM' basm epsi ;
basm = basm 'ET' ('SYME' basm 'DROI' p0 p1) ;
'ELIM' basm epsi ;
* Creation du volume
dx = ri / 2. ;
nz1 = ('ENTIER' ( h1 / dx ))*nn ;
bas = basf 'VOLU' nz1 'TRAN' v1 ;
mbas = basm 'VOLU' nz1 'TRAN' v1 ;
mbas = mbas 'COUL' 'ROUG' ;
plan1 = ax11 'TRAN' nz1 v1 ;
plan4 = ax44 'TRAN' nz1 v1 ;
'ELIM' (bas et plan1 et plan4) epsi ;
mt = bas ;
wall = mbas ;
elim (mt et wall) epsi ;
* Localisation d'une brèche éventuelle au bas de l'enceinte
pjg = 'POIN' basf 'PROC' (0.000 0.000 0.000) ;
jg = ('ELEM' basf 'APPUIE' 'LARGEMENT' pjg) 'COUL' 'VERT' ;
*--------------------------------------------------------------------
* Fin de la définition du maillage
*--------------------------------------------------------------------
*--------------------------------------------------------------------
* Début de l'initialisation de la procédure EXECRXT : table RXT
*--------------------------------------------------------------------
rxt = 'TABLE' ;
*-- Nom du volume fluide
rxt.'vtf' = mt ;
rxt.'epsi' = epsi ;
*-- Definition des murs de l'enceinte : ici un seul mur
*-- en ACIER dont on traitera la thermique dans l'épaisseur
*-- et que l'on initialise a 40oC
*-- On definit d'abort la matériau ACIER avec sa conductivite
*-- thermique LAMBDA (W/m/K) et le produit ro*Cp (J/m3/K)
rxt.'THERMP' = VRAI ;
rxt.'vtp' = wall ;
rxt.'LAMBDA' = 15. ;
rxt.'ROCP' = 3.9E6 ;
rxt.'Tp0' = 60. ;
rxt.'ECHAN' = 10. ;
*-- Conditions initiales dans l'enceinte de test
rxt.'TF0' = 40.0 ;
rxt.'PT0' = 1.00000e5 ;
rxt.'Yvap0' = 0.0023 ;
*-- On positionne une brèche
rxt.'breche' = jg ;
rxt.'diru1' = 0 0 1 ;
*-- On definit un point interne au maillage pour imposer la valeur de
*-- la pression
rxt.'pi' = (0.0 0.0 0.5) ;
*-- On indique que le calcul comporte de la vapeur d'eau
rxt.'VAPEUR' = VRAI ;
*-- On active le recalcul automatique du préconditionnement
*-- toutes les 5 itérations
rxt.'FRPREC' = 5 ;
rxt.'DETMAT' = VRAI ;
*-- Definition du scenario thermohydraulique
rxt.'scenario' = table ;
*-- Conditions a la breche (Obligatoire pour l'instant)
rxt.'scenario'.'t' = prog 0.0 1000.0 ;
rxt.'scenario'.'qeau' = prog 0.050 0.050 ;
rxt.'scenario'.'qair' = prog 0.000 0.000 ;
rxt.'scenario'.'tinj' = prog 150.0 150.0 ;
*-- On impose le pas de temps (s)
rxt.'DT0' = DT0 ;
*-- On impose la viscosite turbulente (m2/s)
rxt.'MODTURB'='NUTURB' ;
rxt.'NUT' = 1.e-2 ;
*-- On lance le calcul sur 20 itérations d'une seconde
rxt.'GRAPH'=GRAPH ;
EXECRXT 2 rxt ;
EXECRXT (nbit - 2) rxt ;
list rxt.TIC.'Tfm' ;
list rxt.TIC.'PT' ;
list rxt.TIC.'Qc' ;
list rxt.TIC.'LMAXU';
ltfm=Prog 40.000 65.434 73.712 81.158 86.800
90.634;
lPT =Prog 1.00000E+05 1.06009E+05 1.21199E+05 1.29592E+05
1.36458E+05 1.43948E+05;
Lqc =Prog 0.0000 0.0000 0.0000 4.41098E-04
3.03436E-03 5.60231E-03;
Lmaxu=Prog 0.0000 0.81320 2.0813 2.6901
2.2350 2.5277;
tic=rxt.'TIC' ;
ERtf=resu( abs (ltfm - tic.'Tfm') )/ 80. ;
ERPT=resu( abs (lPT - tic.'PT' ) ) /1.e5 ;
ERQc=resu( abs (lqc - tic.'Qc' ) ) ;
ERum=resu( abs (Lmaxu - tic.'LMAXU' ) )/2. ;
Mess 'ERtf=' ERtf 'ERPT=' ERPT 'ERQc=' ERQc 'ERum=' ERum ;
Si (ERtf '>' 1.e-4) ; erreur 5 ; Finsi ;
Si (ERPT '>' 1.e-3) ; erreur 5 ; Finsi ;
Si (ERQc '>' 1.e-4) ; erreur 5 ; Finsi ;
Si (ERum '>' 1.e-2) ; erreur 5 ; Finsi ;
Si GRAPH ;
$vtf=rxt.'GEO'.'$vtf' ;
mt =doma $vtf maillage ;
Mpl1=chan 'QUAF' plan1 ;
Mpl4=chan 'QUAF' plan4 ;
ELIM (mt et Mpl1 et Mpl4) epsi ;
$mpl1= mode Mpl1 'NAVIER_STOKES' MACRO ;
$mpl4= mode Mpl4 'NAVIER_STOKES' MACRO ;
plan1= doma $mpl1 maillage ;
plan4= doma $mpl4 maillage ;
plan=plan1 et plan4 ;
paroif = rxt.'GEO'.'paroif';
rho=rxt.'TIC'.'RHO';
rvp=rxt.'TIC'.'RVP';
tf=rxt.'TIC'.'TF';
un=rxt.'TIC'.'UN' ;
un1=redu un plan ;
ung= vect un1 1. ux uy uz jaune ;
trace ung plan ;
opti isov suli ;
trace rho plan 'TITRE'' Rho' ;
trace rvp plan 'TITRE'' Rvp' ;
trace tf plan 'TITRE'' Tf ' ;
trace rho paroif 'TITRE'' Rho' ;
trace rvp paroif 'TITRE'' Rvp' ;
trace Tf paroif 'TITRE'' Tf ' ;
Fcond = tic.'Fcondw';
trace Fcond paroif 'TITRE' ' Flux de condensation Kg / m**2 ' ;
axe = p0 d nz1 (p0 plus v1) ;
axe = chan axe 'QUAF' ;
elim (axe et mt) epsi ;
evr= evol 'CHPO' rho axe ;
dess evr 'TITRE' ' Rho axe ';
evt = evol 'CHPO' Tf axe ;
dess evt 'TITRE' ' Température axe ';
evvp= evol 'CHPO' rvp axe ;
dess evvp 'TITRE' ' rvp axe ';
'FINSI' ;
'FIN' ;
© Cast3M 2003 - All rights reserved.
Disclaimer