Accueil | Cast3M

* fichier :  rayo_abs-2D-2.dgibi
************************************************************************
************************************************************************
 
complet = faux;
*
*************************************************
*  pour calcul complet mettre complet à : vrai;
*************************************************
************************************************************************
*                                                                      *
* Rayonnement thermique en milieu absorbant dans une cavité cylindrique*
* (pas de couplage avec d'autres modes de transfet d'énergie)          *
*                                                                      *
* Comparaison à un calcul analytique                                   *
*                                                                      *
*                                                                      *
* On evalue la puissance perdue par un gaz absorbant de température    *
* uniforme 2273K et de coefficient d'absorption 100./m contenu dans    *
* cavité de rayon 0.01m à la température uniforme de 1273K.            *
*                                                                      *
* Calcul 2D plan                                                       *
* Remarque: le calcul des facteurs de forme n'utilise pas l'option     *
*            convexe 'CONV'                                            *
************************************************************************
 option echo 1 dime 2 elem qua4;
 graph = faux ;
 
*-----------------------------------------------------
* fonction de Bickley-Naylor à l'ordre 3 (normalisée)
* intervenant dans le calcul du facteur de forme
'DEBPROC' BICK3   tau*'FLOTTANT';
 
* intervalle et nombre de points d'intégration
tetat = pi/2.;
nbou = 50 ;
dteta = tetat / nbou;
 
xki1 = 0. ;
repe bcl nbou ;
teta = dteta * (0.5 + (&bcl-1));
sinteta =   'SIN' (teta*90.*2./pi) ;
zz = tau / sinteta ;
zz = -1 * zz ;
dki1 = ('EXP' zz)*dteta;
dki1 = dki1 * sinteta *sinteta ;
xki1 = xki1 + dki1;
fin bcl ;
* normalisation
xki1 = xki1 * 4./pi   ;
*  mess 'bickley xki1: ' xki1 ;
'FINPROC' xki1;
*-----------------------------------------------------
*
*-------------------------------
*   Maillage 
*------------------------------- 
*  
R1 = 1.E-2 ; 
O  =  0. 0. ;
p1 =  0. ( -1. * R1 )  ;  
p2 =  R1  0. ;  
p3 =  0. R1 ;
 
si complet ;
 d = 2.E-4 ;
sinon;
 d = 1.E-3 ;
finsi;
 
p1p2 = cerc p1 O p2 DINI d DFIN d; 
p2p3 = cerc p2 O p3 DINI d DFIN d;
sphe_dro = p1p2 et p2p3 ;
sphe_gau = inve (syme 'DROIT' sphe_dro p1 p3) ;
cavite = sphe_dro et sphe_gau ;
elim (d*1e-5) cavite ;
 
si graph ; 
trac cavite ;
finsi ; 
 
tout = cavite;
 
*------------------------------- 
*  Propriétés physiques 
*------------------------------- 
 
e_wall = 0.5 ;
 abso0 = 100. ;
T_wall = 1273.  ; 
T_gas = 2273. ; 
 
tau1 = 2.*R1*abso0;
 
*-------------------------------
* puissance théorique perdue par les frontieres:
*
* calcul de la solution analytique
*------------------------------- 
 
* calcul du facteur de forme F11 du cylindre sur lui-même
* intervalle et nombre de points d'intégration
 
betat= pi/2. ;
nbeta = 50 ;
dbeta = betat/nbeta;
 
F11 = 0. ;
repe bcl1 nbeta ;
beta = dbeta * (0.5 + (&bcl1-1));
cosbeta = 'COS' (beta*90.*2./pi) ;
tau2 = tau1 * cosbeta;
Ki1 = BICK3   tau2 ;
*mess ' Ki3 : ' Ki1;
DF = Ki1 *cosbeta*dbeta;
F11 = F11 + DF ;
fin bcl1;
 
* facteur de forme
 
mess 'F11: ' F11 ;
mess '1-F11: ' (1-F11) ;
 
* puissance échangée cylindre-milieu absorbant
 
b0 = 5.67e-8;
EMG1 = b0 * ( T_gas ** 4. ) ;
EMS1 = b0 * ( T_wall ** 4. ) ;
AIRE1 = 2. *pi * R1 ;
abso1 = 1. - tau1 ;
*denom0 = ( 1. / e_wall ) + ( 1. / abso1 ) - 1. ;
res0 = -1. * AIRE1 * ( EMG1 - EMS1 ) ;
coef1 = 1. - F11 ;
coef2 = 1. - ((1. - e_wall)*F11) ;
coef = e_wall * coef1 /coef2 ;
res0 = coef * res0 ;
mess 'res0: ' res0 ;
 
*------------------------------- 
*  Modèle de rayonnement 
*------------------------------- 
 
mrt = mode cavite thermique   rayonnement 'CAVITE' ;
e   = mate mrt 'EMIS' e_wall  CABS (-1 *abso0) 'TABS' T_gas ; 
 
*------------------------------- 
*  Facteurs de forme et matrice de rayonnement
*------------------------------- 
opti impi 1; comm pour récupérer le bilan ;
 fft = ffor mrt e    ;
opti impi 0;
 
 chamr = raye mrt  fft    e ;
 
*------------------------------- 
* gaz absorbant : calcul du terme R*Tg4
*------------------------------- 
 
tg = manu chpo tout 1 'T'  T_gas natu 'DIFFUS';
tg_cavi   = redu tg cavite ;
tge_cavi   =  chan 'CHAM' tg_cavi   mrt 'GRAVITE' ;
crg= rayn mrt  chamr tge_cavi    ;
fg =  crg*tg_cavi ;
 
*-------------------------------  
* paroi : calcul du terme R*Tw4
*------------------------------- 
 
tp = manu chpo tout 1 'T' T_wall natu 'DIFFUS';
t_cavi   = redu tp cavite ;
te_cavi   =  chan 'CHAM' t_cavi   mrt 'GRAVITE' ;
cr = rayn mrt  chamr te_cavi    ;
fw =  cr *t_cavi ;
 
*-------------------------------
* puissance perdue par les frontieres:
*------------------------------- 
 
puis_n =  fw - fg ;
t1 = manu chpo tout 1 'T' 1.0    natu 'DIFFUS' ;
lm1 = MOTS 'Q' ;  
lm2 = MOTS 'T' ;
puis_t = xty puis_n t1 lm1 lm2 ;
 
 
*------------------------------- 
*  Erreur methode 1
*------------------------------- 
 
si ( ( 'ABS' res0 ) '>' 1.E-5 ) ; 
  err1 = ( 'ABS' ( res0 - puis_t ) ) / ( 'ABS' res0 ) ;
  err1 = err1 * 100. ;
sinon ; 
  err1 = 0. ;
finsi ; 
'MESS' ' Puissance absorbée calculée methode 1  = ' puis_t ;
'MESS' ' Puissance absorbee théorique = ' res0 ;
'MESS' ' Erreur methode 1 en pourcentage = ' err1 ' % ' ;
 
 
*------------------------------- 
* calcul du flux rayonné au moyen de l'évaluation de la
* température de rayonnement (option 2 de RAYE)
*------------------------------- 
 
       trad = raye mrt  fft e  te_cavi 1e-7  ;
*      mess 'trad ' (mini trad) (maxi trad);
 
* calcul du coefficient d'echange
 
       hrad = HRCAV  mrt e te_cavi trad ;
*      mess 'hrad ' (mini hrad) (maxi hrad);
 
       trad_n1= chan 'CHPO' mrt trad       ;
       trad_n = nomc trad_n1 'T' 'NATU' 'DIFFUS';
 
* pour la condition de convection
       cr = cond mrt hrad;
       f  = conv mrt hrad trad_n;
 
*      flux rayonne
 
       fray  =  (cr * tp)- f;
       fray_tot = maxi (resu fray);
 
       mess ' flux par methode 2 ' fray_tot ;
 
*------------------------------- 
*  Erreur methode 2
*------------------------------- 
 
si ( ( 'ABS' res0 ) '>' 1.E-5 ) ; 
  err2 = ( 'ABS' ( res0 - fray_tot) ) / ( 'ABS' res0 ) ;
  err2 = err2 * 100. ;
sinon ; 
  err2 = 0. ;
finsi ; 
'MESS' ' Puissance absorbée calculée methode 2  = ' fray_tot ;
'MESS' ' Puissance absorbee théorique = ' res0 ;
'MESS' ' Erreur methode 2 en pourcentage = ' err2 ' % ' ;
 
 
si (( err1 '<' .2 ) et ( err2 '<' .2 ));
 'ERRE' 0 ;
sinon ; 
 'ERRE' 5 ;
finsi ; 
 
fin ;