Télécharger la plaquette Cast3M


fr en

Télécharger Marangoni3.dgibi

Retour à la liste

Numérotation des lignes : 

  1. * fichier :  Marangoni3.dgibi
  2. *                      _____________________
  3. *
  4. *        CAS TEST D'UN ECOULEMENT SOUS UN FLUX D'ENERGIE CONCENTRE
  5. *             UTILISATION DE LA LOI DE CHAN, CHEN, MAZUMDER
  6. *           Journal of Heat Transfert, vol 110, 1988
  7. *                      _____________________
  8.  
  9. ***************************************************************
  10. ** DESCRIPTION DU TEST :
  11. **
  12. ** Une poche de Cerium liquide est soumise au bombardement d'un flux de
  13. ** chaleur surfacique ayant une loi parabolique. La loi de CHAN décrit
  14. ** ce type d'écoulement. Cette loi  donne le profil de température du
  15. ** liquide ainsi que le profil des vitesses a la surface de la poche
  16. ** liquide.
  17. **
  18. ** Lors de ce test, nous imposons donc un flux surfacique sur une poche
  19. ** de metal liquide. Nous comparons les profils obtenus de température
  20. ** et de vitesse a la surface du liquide avec les valeurs théoriques 
  21. ** données par la loi de CHAN.
  22. ***************************************************************
  23.  
  24. option DIME 2 ELEM qua8 MODE axis; MOAXI = VRAI ;
  25. *option DIME 2 ELEM qua8 MODE plan; MOAXI = FAUX ;
  26.  
  27. option trac psc ;
  28.  
  29. *                 _________________________________
  30. *
  31. *                  différents parametres du calcul
  32. *                 _________________________________
  33.  
  34. COMPLET=FAUX ;
  35. GRAPH  =FAUX ;
  36. * Nombre de mailles
  37. NMAIL = 20 ;
  38.  
  39. * nombre maximal d'itérations
  40. 'SI' COMPLET ;NITER = 200 ;
  41. 'SINON' ;NITER = 4 ;
  42. 'FINSI' ;
  43. * coefficient de relaxation pour le champ de température
  44. ALFAT = 0.2 ;
  45.  
  46. * coefficient de relaxation pour le champ de vitesse
  47. ALFAU = 0.5 ;
  48.  
  49. * précision requise
  50. PREC = 5.e-3 ;
  51.  
  52. * pseudo pas de temps pour la qdm
  53. dt = .01 ;
  54.  
  55. *   _______________________________________________________
  56. *
  57. *    géométrie du creuset : tronc de cone (dimension en m)
  58. *   _______________________________________________________
  59. *
  60.  
  61. * rayon du bas du creuset
  62. rbas  = 0.0399 ;
  63.  
  64. * rayon du haut du creuset
  65. rhaut = 0.04 ;
  66.  
  67. * hauteur de liquide
  68. hliq  = 0.020 ;
  69.  
  70. * épaisseur de la couche limite
  71. eclm = 1.e-4 ;
  72.  
  73. * on ne maille pas l'axe, à cause des termes en 1/r et 1/r2
  74. ref = 0.    ;
  75.  
  76. * rayon moyen du faisceau électonique
  77. rbe   = 0.010 ;
  78.  
  79. *  __________________________________________________________________
  80. *
  81. *   Caractéristiques du faisceau d'électrons à répartition d'énergie
  82. *    gaussienne : puissance PCANON et "rayon" RFOCAL
  83. *  __________________________________________________________________
  84.  
  85. PCANON = 1.e3 ;
  86. RFOCAL = rbe ;
  87.  
  88.  
  89. *                 ______________________
  90. *
  91. *                  constantes physiques
  92. *                 ______________________
  93.  
  94. * constante des gaz parfaits
  95. RG = 8.32 ;
  96.  
  97. * constante de STEPHAN-BOLTZMAN
  98. SIGMA = 5.67E-8 ;
  99.  
  100. * accélération de la pesanteur
  101. g = -9.81 ;
  102.  
  103. *
  104. *         __________________________________________________
  105. *
  106. *          propriétés physiques du cérium évaluées à 2500 K
  107. *         __________________________________________________
  108.  
  109. * viscosité dynamique (kg/m/s)
  110. mu = 6.5e-3     ;
  111.  
  112. * masse volumique (kg/m3)
  113. rho0 = 6310. ;
  114.  
  115. * conductivité thermique  (W/m/K)
  116. lambda = 25. ;
  117.  
  118. * chaleur spécifique (J/kg/K)
  119. Cp = 270 ;
  120.  
  121. * température de fusion (K)
  122. tfus = 1070. ;
  123.  
  124. * cooefficient de dilatation volumique (1/K)
  125. beta = 4.2e-5 ;
  126.  
  127. * tension de surface (N/m)
  128. sig0 = 0.34 ;
  129.  
  130. * 1/sig0 dsig0/dT (1/K)
  131. gam0 = 5.88e-4 ;
  132.  
  133. * émissivité
  134. emis = 0.34 ;
  135.  
  136. * coefficient d'absorption des électrons (% d'énergie absorbée)
  137. CABSOR = 1. ;
  138.  
  139. * masse molaire (kg/mole)
  140. MM = 0.14012 ;
  141.  
  142. * chaleur latente de vaporisation (J/kg)
  143. LAT = 2.95e6 ;
  144.  
  145. * coefficients pour le calcul de la pression de vapeur saturante
  146. PA = 6.023 ;
  147. PB = -21278. ;
  148. PC = -0.1127 ;
  149.  
  150. * coefficient de recondensation de la vapeur
  151. RET = 0. ;
  152.  
  153. * viscosité cinématique (m2/s)
  154. nu = mu / rho0 ;
  155.  
  156. * diffusivité thermique (m2/s)
  157. khi = lambda / ( rho0 * Cp ) ;
  158.  
  159. * effet MARANGONI
  160. lev =  1.  * gam0 * sig0 / rho0 ;
  161.  
  162. *                     ______________________
  163. *
  164. *                      CREATION DU MAILLAGE
  165. *                     ______________________
  166.  
  167. p0  = ref    0.          ;
  168. p1  = rhaut  0.          ;
  169. p2  = rbas   (-1. * hliq) ;
  170. p3  = ref    (-1. * hliq) ;
  171.  
  172. p01 = rbe    0.          ;
  173. p32 = rbe    (-1. * hliq) ;
  174.  
  175. p03 = ref    (-10. * eclm) ;
  176. p12 = rhaut  (-10. * eclm) ;
  177.  
  178. d0 = rbe / NMAIL ;
  179. d1 = rbe * 4 / NMAIL ;
  180. d2 = rhaut / 2 / NMAIL ;
  181. d3 = rhaut * 6 / NMAIL ;
  182. surli0 = (p01 droi dini d1 dfin d0 p0 ) coul roug ;
  183. surli1 = (p1  droi dini d3 dfin d2 p01) coul roug ;
  184. surli = surli1 et surli0 ;
  185.  
  186. d0 = eclm * 20 / NMAIL ;
  187. d1 = eclm * 60 / NMAIL ;
  188. d2 = hliq / NMAIL ;
  189. d3 = hliq * 8 / NMAIL ;
  190. smm0 = (p0  droi dini d0 dfin d1 p03) coul jaun ;
  191. smm1 = (p03 droi dini d2 dfin d3 p3 ) coul jaun ;
  192. smm = smm1 et smm0 ;
  193.  
  194. scon10 = (inve (smm0 proj 'DIRE' (1. 0.) 'DROI' p1 p2) ) coul vert ;
  195. scon11 = (inve (smm1 proj 'DIRE' (1. 0.) 'DROI' p1 p2) ) coul vert ;
  196. scon1 = scon11 et scon10 ;
  197.  
  198. dz0 = (smm0 proj 'DIRE' (1. 0.) 'DROI' p01 p32) coul bleu ;
  199. dz1 = (smm1 proj 'DIRE' (1. 0.) 'DROI' p01 p32) coul bleu ;
  200. dz = dz0 et dz1 ;
  201.  
  202. zc = (-1. * (rhaut - rbe ) * hliq ) / ( rhaut - rbas) ;
  203.  
  204. scon21= (inve (surli1 proj 'CONI' (rbe zc) 'DROI' p2 p3) ) coul vert ;
  205. scon20= (inve (surli0 proj 'DIRE' (0. -1.) 'DROI' p2 p3) ) coul vert ;
  206. scon2 = scon20 et scon21 ;
  207.  
  208. dr1 = (surli1 proj 'CONI' (rbe zc) 'DROI' p12 p03) coul rose ;
  209. dr0 = (surli0 proj 'DIRE' (0. -1.) 'DROI' p12 p03) coul rose ;
  210. dr = dr1 et dr0 ;
  211.  
  212. scon = scon1 et scon2 ;
  213.  
  214. elim 1.e-6 ( surli1 et surli0 et scon10 et scon11 et scon20
  215.            et scon21 et dr0 et dr1 et dz0 et dz1 et smm0 et smm1 ) ;
  216.  
  217. meta1 = daller surli1 dz0 dr1 scon10 ;
  218. meta2 = daller (inve dr1) dz1 scon21 scon11 ;
  219. meta3 = daller (inve dz0) surli0 smm0 dr0 ;
  220. meta4 = daller (inve dz1) (inve dr0) smm1 scon20 ;
  221.  
  222. meta = meta1 et meta2 et meta3 et meta4 ;
  223.  
  224.  
  225. oublier meta1 ;
  226. oublier meta2 ;
  227. oublier meta3 ;
  228. oublier meta4 ;
  229. 'SI' GRAPH ;trac meta ; 'FINSI' ;
  230.  
  231. *
  232. ** On réoriente les éléments.
  233. *
  234.  
  235. meta = 'ORIENTER' meta ;
  236. tass meta ;
  237. nbel1 = nbel meta ; list nbel1 ;
  238.  
  239. * transformation des éléments en qua9
  240.  
  241. £META  = CHANGER META  MACRO ;
  242. £SURLI = CHANGER SURLI MACRO ;
  243. £SURLI0 = CHANGER SURLI0 MACRO ;
  244. £SCON  = CHANGER SCON  MACRO ;
  245. £SCON2 = CHANGER SCON2 MACRO ;
  246. £SMM   = CHANGER SMM   MACRO ;
  247.  
  248. ELIM 1.e-6 (£META ET £SURLI ET SURLI0 ET £SCON ET SCON2 ET £SMM) ;
  249.  
  250. * formulation du modèle NAVIER_STOKES
  251.  
  252. $META  = MODE £META  'NAVIER_STOKES' MACRO ;
  253. $SURLI = MODE £SURLI 'NAVIER_STOKES' MACRO ;
  254. $SURLI0 = MODE £SURLI0 'NAVIER_STOKES' MACRO ;
  255. $SCON  = MODE £SCON  'NAVIER_STOKES' MACRO ;
  256. $SCON2 = MODE £SCON2 'NAVIER_STOKES' MACRO ;
  257. $SMM   = MODE £SMM   'NAVIER_STOKES' MACRO ;
  258.  
  259.  
  260. mmeta = doma $META  maillage ;
  261. ssurli = doma $SURLI maillage ;
  262. ssurli0 = doma $SURLI0 maillage ;
  263.  
  264. *       _______________________________________________
  265. *
  266. *        SYSTEME D'EQUATIONS ET CONDITIONS AUX LIMITES
  267. *       _______________________________________________
  268.  
  269.  
  270.  
  271.  
  272. * equations de qdm sur ur et uz + div V = 0
  273.  
  274. RV = EQEX $META
  275. OPTI 'EF'  'IMPL' CENTREP1  'SUPG'
  276. ZONE $META OPER KBBT 1 (-1.)                INCO 'UN' 'PN'
  277. ZONE $META OPER NS 1. 'UN' nu               INCO 'UN'
  278. ZONE $META OPER DFDT 1.  'UN'  dt 'UN'  nu  INCO 'UN'
  279. OPTI 'CENTREE'
  280. ZONE $SURLI OPER TOIMP 'tau'                INCO 'UN'
  281. ;
  282.  
  283.  
  284. * equation scalaire sur la température + flux imposé sur surli
  285.  
  286. RV = EQEX RV
  287. OPTI 'EF'  'IMPL' 'SUPG'
  288. ZONE $META  OPER TSCA 1. 'UN' khi 0.        INCO 'TN'
  289. ZONE $SURLI OPER FIMP 'F1E'               INCO 'TN'
  290. ;
  291.  
  292. * conditions aux limites
  293.  
  294.  
  295. RV = EQEX RV
  296. CLIM UN UIMP  £SMM     0.
  297.      UN UIMP £SCON2    0.
  298.      UN VIMP £SURLI    0.
  299.      TN TIMP £SCON2    tfus
  300. ;
  301.  
  302. *                      ________________
  303. *
  304. *                       INITIALISATION
  305. *                      ________________
  306.  
  307. RV.INCO = TABLE INCO ;
  308.  
  309. RV.INCO.'UN'   = KCHT $META VECT SOMMET (1.e-8 1.e-8) ;
  310. RV.INCO.'TN'   = KCHT $META SCAL SOMMET tfus ;
  311. RV.INCO.'PN'   = KCHT $META SCAL CENTREP1 0. ;
  312. RV.INCO.'tau'  = KCHT $SURLI VECT CENTRE (0. 0. ) ;
  313.  
  314.  
  315. * terme f1e  flux imposé sur la surface libre
  316.  
  317. I0 = PCANON /(PI * (RFOCAL**2.)) ;
  318. XK = ssurli coor 1 ;
  319. XK = KCHT $SURLI SCAL SOMMET XK ;
  320. FBA = 'KOPS' I0 * ('KOPS' 1. - ( 'KOPS' ('KOPS' XK '/' rbe) ** 2.) );
  321. mask1 = XK MASQUE 'EGINF' rbe ;
  322. mask1 = 'KCHT' $SURLI SCAL SOMMET mask1 ;
  323. FBA = 'KOPS' FBA * mask1 ;
  324. FBA = KCHT $SURLI SCAL SOMMET FBA ;
  325. F1 = 'KOPS' FBA '/' ( rho0 * Cp ) ;
  326. F1E = NOEL $SURLI F1 ;
  327. RV.INCO.'F1E' = KCHT $SURLI SCAL CENTRE F1E ;
  328.  
  329. *                     __________________
  330. *
  331. *                      BOUCLE ITERATIVE
  332. *                     __________________
  333.  
  334. REPETER BLOCKI NITER ;
  335.  
  336. *
  337. ** calcul des termes de couplages
  338. *
  339.  
  340. MESS 'ITERATION INTERNE N°' &BLOCKI ;
  341.  
  342. UN0 = COPI RV.INCO.'UN' ;
  343. TN0 = COPI RV.INCO.'TN' ;
  344.  
  345.  
  346.  
  347. * terme tau
  348. * grad T < 0 ,lev > 0 => Marangoni < 0  donc opposé à la qdm
  349. *
  350. TK = REDU TN0 ssurli ;
  351. TK = KCHT $SURLI SCAL SOMMET TK ;
  352. GTK = KOPS GRAD TK $SURLI ;
  353. mask1 = GTK MASQUE 'EGINF' 0. ;
  354. mask1 = 'KCHT' $SURLI VECT CENTRE mask1 ;
  355. GTK = 'KOPS' GTK * mask1 ;
  356. MARK = KOPS lev '*' GTK ;
  357. RV.INCO.'tau' = KCHT $SURLI VECT CENTRE MARK ;
  358. MARKN = ELNO $SURLI MARK ;
  359.  
  360.  
  361. *
  362. * CALCUL D'UNE ITERATION
  363. *
  364. RV.METHINV.FCPRECT = 1 ;
  365. RV.METHINV.FCPRECI = 1 ;
  366. EXEC RV ;
  367.  
  368.  
  369. *
  370. *                      TESTS DE CONVERGENCE
  371. *
  372. UNR = EXCO UX RV.INCO.'UN' ;
  373. UNZ = EXCO UY RV.INCO.'UN' ;
  374. URMAX = MAXI UNR ;
  375. UZMAX = MAXI UNZ ;
  376. TMAX = MAXI RV.INCO.'TN' ;
  377. MESS 'TMAX = ' TMAX ;
  378. MESS 'URMAX = ' URMAX ;
  379. MESS 'UZMAX = ' UZMAX ;
  380.  
  381.  
  382.  
  383. UR0 = UN0 EXCO 'UX' scal ;
  384. UR0 = KCHT $META SCAL SOMMET UR0 ;
  385. UNR = KCHT $META SCAL SOMMET UNR ;
  386. ERUNR = KOPS ( KOPS UNR '-' UR0 ) '/' ( (MAXI (ABS UNR) ) + 1.e-13) ;
  387. ERRUNR = MAXI ( ABS ERUNR ) ;
  388. MESS 'ERRUNR = ' ERRUNR ;
  389.  
  390. UZ0 = EXCO UY UN0 ;
  391. UZ0 = KCHT $META SCAL SOMMET UZ0 ;
  392. UNZ = KCHT $META SCAL SOMMET UNZ ;
  393. ERUNZ = KOPS ( KOPS UNZ '-' UZ0 ) '/' ( (MAXI (ABS UNZ) ) + 1.e-13) ;
  394. ERRUNZ = MAXI ( ABS ERUNZ ) ;
  395. MESS 'ERRUNZ = ' ERRUNZ ;
  396.  
  397.  
  398. TN = RV.INCO.'TN' ;
  399. ERTN = KOPS ( KOPS TN '-' TN0 ) '/' ( (MAXI (ABS TN) ) + 1.e-13) ;
  400. ERRTN = MAXI ( ABS ERTN ) ;
  401. MESS 'ERRTN = ' ERRTN ;
  402.  
  403. ERR = ERRUNR ; VAR = MOTS 'UNR' ;
  404.  
  405. SI ( ERRUNZ > ERR ) ;
  406.   ERR = ERRUNZ ; VAR = MOTS 'UNZ' ;
  407. FINSI ;
  408.  
  409.  
  410. SI ( ERRTN > ERR ) ;
  411.   ERR = ERRTN ; VAR = MOTS 'TN' ;
  412. FINSI ;
  413.  
  414. VAR1 = EXTR 1 VAR ;
  415. MESS 'erreur relative maximale sur ' VAR1 ERR ;
  416.  
  417. * RELAXATION
  418.  
  419. RV.INCO.'TN' = KOPS ( KOPS ALFAT '*' RV.INCO.'TN')
  420.                     '+' ( KOPS (1. - ALFAT) '*' TN0 ) ;
  421.  
  422. RV.INCO.'UN' = KOPS ( KOPS ALFAU '*' RV.INCO.'UN')
  423.                     '+' ( KOPS (1. - ALFAU) '*' UN0 ) ;
  424.  
  425.  
  426. SI (ERR &lt;EG PREC) ;
  427. * convergence des itérations internes
  428.   MESS 'CONVERGENCE EN ' &BLOCKI 'ITERATIONS' ;
  429.   QUITTER BLOCKI ;
  430. FINSI ;
  431.  
  432. SI (&BLOCKI EGA NITER ) ;
  433.   MESS 'PAS DE CONVERGENCE EN ' NITER 'ITERATIONS' ;
  434. FINSI ;
  435.  
  436.  
  437.  
  438. ******************************
  439. ** FIN DE LA BOUCLE ITERATIVE
  440. *******************************
  441. FIN BLOCKI ;
  442.  
  443. ******************************
  444. ** COURBE THEORIQUE
  445. ******************************
  446.  
  447. 'SI' ( MOAXI ) ;
  448. bb = 0.70185 ;
  449. cc = 0.35093 ;
  450. f1 = 0.89872 ;
  451. SINON ;
  452. bb = 1. ;
  453. cc = 0.33333 ;
  454. f1 = 1. ;
  455. 'FINSI' ;
  456.  
  457. q0 = I0 ;
  458. q1 = q0 '/' rbe '/' rbe ;
  459. XK = ssurli coor 1 ;
  460.  
  461. ** Courbe theorique de temperature sur la surface selon CHAN
  462. DTmax = bb * (cc ** -.25) * q0 * (q1 ** -.25) 
  463.         * ((2 * (rho0 ** 2) * (Cp ** 3) * lev * mu) ** -.25) ;
  464. 'LISTE' DTmax ;
  465. Ttheo = 'KOPS' DTmax * ('KOPS' 1 '-' ('KOPS' (cc / bb / rbe / rbe)
  466.          * ('KOPS' XK ** 2 ) ) ) ;
  467. Ttheo = 'KOPS' Ttheo '+' Tfus ;
  468. Ttheo = KCHT $SURLI SCAL SOMMET Ttheo ; 
  469.  
  470. ** Courbe theorique de la vitesse radiale sur la surface selon CHAN
  471. URtheo = ((cc * q1 * 2 * lev / mu / Cp) ** .5) * f1 ;
  472. URtheo = 'KOPS' URtheo * XK ;
  473. URtheo = KCHT $SURLI SCAL SOMMET URtheo ; 
  474. evURtheo = evol roug chpo URtheo (inve ssurli ) ;
  475.  
  476. ******************
  477. ** IMPRESSIONS
  478. ******************
  479. titre 'Temperature le long de surli0 (W/m2)' ;
  480. TN = RV.INCO.'TN' ;
  481. evtemp = evol noir chpo TN (inve ssurli0 ) ;
  482. evTtheo = evol roug chpo Ttheo (inve ssurli0 ) ;
  483. 'SI' GRAPH ;dess (evtemp 'ET' evTtheo) 'LEGE' ; 'FINSI' ;
  484.  
  485. titre 'vitesses (ampl= 0.01) (m/s)' ;
  486. UN = RV.INCO.'UN' ;
  487. VIT1 = VECT UN 0.01  UX UY ROUG ;
  488. 'SI' GRAPH ;trac VIT1 mmeta ; 'FINSI' ;
  489.  
  490. titre 'vitesse radiale le long de surli (m/s) (en rouge la theorie)' ;
  491. evux = evol chpo UN UX (inve ssurli ) ;
  492. 'SI' GRAPH ;dess (evux 'ET' evURtheo) 'LEGE' ;  'FINSI' ;
  493.  
  494.  
  495. titre 'pression (Pa)' ;
  496. PNN = ELNO $META (RV.'INCO'.'PN') CENTREP1 ;
  497. 'SI' GRAPH ;trace PNN mmeta ; 'FINSI' ;
  498.  
  499. ***************************************
  500. **** TEST DE FONCTIONNEMENT ***********
  501. ***************************************
  502. TN = RV.INCO.'TN' ;
  503. TNmax = 'EXTR' TN 'SCAL' (p0) ;
  504. TNtmax = 'EXTR' Ttheo 'SCAL' (p0) ;
  505. 'LISTE' TNmax ;
  506. 'LISTE' TNtmax ;
  507. ERtheo = ( TNmax - TNtmax ) / ( TNmax - Tfus) ;
  508. ERRtheo = ABS ERtheo ;
  509. MESS 'erreur de la temperature max par rapport à la théorie =' ERRtheo ;
  510.  
  511. 'SI' COMPLET ;
  512.  'SI' (ERRtheo < .05) ;
  513.     'ERREUR' 0 ;
  514.  'SINON' ;
  515.     'ERREUR' 5 ;
  516.  'FINSI' ;
  517. 'SINON' ;
  518.  'SI' (ERRtheo < .4 ) ;
  519.     'ERREUR' 0 ;
  520.  'SINON' ;
  521.     'ERREUR' 5 ;
  522.  'FINSI' ;
  523. 'FINSI' ;     
  524.  
  525. ***************************************
  526.  
  527.  
  528. 'FIN' ;
  529.  
  530.  
  531.  
  532.  
  533.  
  534.  
  535.  
  536.  
  537.  
  538.  
  539.  
  540.  
  541.  
  542.  
  543.  

© Cast3M 2003 - Tous droits réservés.
Mentions légales