Télécharger smithhutton_impl.dgibi

Retour à la liste

Numérotation des lignes :

  1. * fichier : smithhutton_impl.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4. *
  5. *-------------------------------------------------------------------------------
  6. * Convection/Diffusion : CAS SMITH ET HUTTON
  7. * REFERENCE : NUMERICAL HEAT TRANSFER, VOL.5, p.439, 1982
  8. *-------------------------------------------------------------------------------
  9. * Les faces latérales et supérieure d'une boite rectangulaire sont
  10. * imperméables. Le fluide rentre et sort de la boite par la face
  11. * inférieure. Il rentre par la moitié gauche de la face inférieure
  12. * et sort par la moitié droite de cette même face.
  13. *
  14. * Le champ de vitesse est connu et donné par :
  15. * v(x,y) = 2y(1-x2) i - 2x(1-y2) j
  16. * la taille de la boite étant [-1,1]x[0,1]
  17. *
  18. * L'écoulement transporte un scalaire passif qui est injecté en entrée
  19. * suivant le profil :
  20. * c(x,0) = 1 + tanh(10(2x+1))
  21. * avec x variant de -1 à 0 (entrée du domaine).
  22. * Sur les frontières imperméables du domaine, la concentration est
  23. * constante et égale à 1-tanh(10).
  24. *
  25. * On cherche la solution stationnaire du transport par diffusion et
  26. * convection du champ scalaire passif. On compare la solution en sortie
  27. * avec la solution de référence. Plusieurs solutions sont calculées
  28. * suivant le Peclet (rapport entre la convection et la diffusion).
  29. *-------------------------------------------------------------------------------
  30. * Résolution (implicite) de l'équation stationnaire
  31. *-------------------------------------------------------------------------------
  32. *
  33. 'OPTI' 'DIME' 2 'ISOV' 'SULI' ;
  34. *
  35. *
  36. *========================================================================
  37. * USER DATA - USER DATA - USER DATA - USER DATA - USER DATA - USER DATA
  38. *========================================================================
  39. *
  40. *
  41. * Pe : Peclet (convection sur diffusion)
  42. * Les cas traités dans la référence sont 10, 100, 500, 1000 et 1000000
  43. Pe = 1000000 ;
  44. *
  45. * KSUPG : Option de décentrement (CENTREE/SUPG/SUPGDC)
  46. * GRAPH : Booleen pour l'affichage des tracés à l'issue du calcul
  47. * NX : Nombre de maille suivant x
  48. * NY : Nombre de maille suivant y
  49. * tolera : Tolérance pour la comparaison avec la solution de référence
  50. KSUPG = 'CENTREE' ;
  51. GRAPH = faux ;
  52. 'OPTI' 'ELEM' 'QUA4' ;
  53. NY = 20 ;
  54. NX = 2 * NY ;
  55. tolera = 0.0002 ;
  56. *
  57. *========================================================================
  58. * USER DATA - USER DATA - USER DATA - USER DATA - USER DATA - USER DATA
  59. *========================================================================
  60. *
  61. *
  62. *- MAILLAGE
  63. *
  64. *
  65. * Points
  66. A1 = -1.0 0.0 ;
  67. A2 = 1.0 0.0 ;
  68. A3 = 1.0 1.0 ;
  69. A4 = -1.0 1.0 ;
  70. A0 = 0.0 0.0 ;
  71. *
  72. * Lignes
  73. LIN = A1 'DROI' NY A0 ;
  74. LOUT = A0 'DROI' NY A2 ;
  75. FBAS = LIN 'ET' LOUT ;
  76. FDRO = A2 'DROI' NY A3 ;
  77. FHAU = A3 'DROI' NX A4 ;
  78. FGAU = A4 'DROI' NY A1 ;
  79. LIMP = FDRO 'ET' FHAU 'ET' FGAU ;
  80. *
  81. * Maillage
  82. DOMTOT = 'DALL' FBAS FDRO FHAU FGAU 'PLAN' ;
  83. CNT1 = 'CONT' DOMTOT ;
  84. *
  85. * Modèles et sous-modèles
  86. DOM2 = 'CHAN' 'QUAF' DOMTOT ;
  87. LIN2 = 'CHAN' 'QUAF' LIN ;
  88. LIMP2 = 'CHAN' 'QUAF' LIMP ;
  89. $DOMTOT = 'MODE' DOM2 'NAVIER_STOKES' 'LINE' ;
  90. $LIN = 'MODE' LIN2 'NAVIER_STOKES' 'LINE' ;
  91. $LIMP = 'MODE' LIMP2 'NAVIER_STOKES' 'LINE' ;
  92. *
  93. * Récupération des maillages et fusion des supports
  94. DOMTOT = 'DOMA' $DOMTOT 'MAILLAGE' ;
  95. LIN = 'DOMA' $LIN 'MAILLAGE' ;
  96. LIMP = 'DOMA' $LIMP 'MAILLAGE' ;
  97. FDOMTOT = 'DOMA' $DOMTOT 'FACE' ;
  98. CLIN = 'DOMA' $LIN 'CENTRE' ;
  99. CLIMP = 'DOMA' $LIMP 'CENTRE' ;
  100. 'ELIM' FDOMTOT (CLIN 'ET' CLIMP) 1.D-3 ;
  101. *
  102. * Champ de vitesse
  103. XX YY = 'COOR' DOMTOT ;
  104. VXSH = (2.*YY)*(1.0-(XX*XX)) ;
  105. VYSH = (-2.*XX)*(1.0-(YY*YY)) ;
  106. VX0 = 'NOMC' 'UX' VXSH ;
  107. VY0 = 'NOMC' 'UY' VYSH ;
  108. VX = 'KCHT' $DOMTOT 'SCAL' 'SOMMET' 'COMP' 'UX' VX0 ;
  109. VY = 'KCHT' $DOMTOT 'SCAL' 'SOMMET' 'COMP' 'UY' VY0 ;
  110. CHVIT = 'KCHT' $DOMTOT 'VECT' 'SOMMET' 'COMP' 'UX' 'UY' (VX 'ET' VY) ;
  111. *
  112. * Diffusion (1/Pe)
  113. DIF = 1.0 / ('FLOT' Pe) ;
  114. *
  115. * Profil de concentration à l'entrée
  116. XBAS = 'COOR' 1 LIN ;
  117. TOTO = 2.0*XBAS + 1.0 * 10. ;
  118. SOLUTION = 'TANH' TOTO ;
  119. SOLUTION = 1.0 + SOLUTION ;
  120. CHP1 = 'KCHT' $LIN 'SCAL' 'SOMMET' 0. SOLUTION ;
  121. CHP1 = 'NOMC' 'CN' CHP1 ;
  122. *
  123. * Conditions aux limites en concentration sur les frontières imperméables
  124. C1 = (1.0 - (TANH 10.0)) ;
  125. *
  126. * Description du problème de transport
  127. RV1 = 'EQEX' $DOMTOT 'ITMA' 1 'ALFA' 0.7
  128. 'OPTI' 'EF' 'IMPL' KSUPG
  129. * 'ZONE' $DOMTOT 'OPER' 'TSCAL' DIF 'VITESSE' 0. 'INCO' 'CN'
  130. 'ZONE' $DOMTOT 'OPER' 'LAPN' DIF 'INCO' 'CN'
  131. 'ZONE' $DOMTOT 'OPER' 'KONV' 1. 'VITESSE' DIF 'INCO' 'CN'
  132. ;
  133. *
  134. * Description des conditions aux limites
  135. RV1 = 'EQEX' RV1
  136. 'CLIM' 'CN' 'TIMP' LIN CHP1
  137. 'CLIM' 'CN' 'TIMP' LIMP C1 ;
  138. *
  139. * Description des conditions initiales
  140. RV1 . 'INCO' = TABLE 'INCO' ;
  141. RV1 . 'INCO' . 'CN' = 'KCHT' $DOMTOT 'SCAL' 'SOMMET' 0. ;
  142. *
  143. * Autres data
  144. RV1 . 'INCO' . 'VITESSE' = CHVIT ;
  145. RV1 . 'INCO' . 'CN2' = 'KCHT' $DOMTOT 'SCAL' 'SOMMET' 0. ;
  146. RV1 . 'INCO' . 'IT' = 'PROG' ;
  147. RV1 . 'INCO' . 'TI' = 'PROG' ;
  148. RV1 . 'INCO' . 'ER' = 'PROG' ;
  149. *
  150. *
  151. *- CALCUL
  152. *
  153. *
  154. EXEC RV1 ;
  155. *
  156. *
  157. *- ANALYSE DES RESULTATS
  158. *
  159. *
  160. EVOL1 = 'EVOL' 'CHPO' (RV1 . 'INCO' . 'CN') 'SCAL' FBAS ;
  161. EVOL2 = 'EVOL' 'CHPO' XX 'SCAL' FBAS ;
  162. LIX = 'EXTR' EVOL2 'ORDO' ;
  163. LIU = 'EXTR' EVOL1 'ORDO' ;
  164. EVOL3 = 'EVOL' 'MANU' 'X' LIX 'C(X,0)' LIU ;
  165. *
  166. * Solutions de référence
  167. LIXT = PROG 0. 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 ;
  168. SI ( Pe 'EGA' 10 ) ;
  169. LIUT = PROG 1.989 1.402 1.146 0.946 0.775 0.621 0.480 0.349
  170. 0.227 0.111 0.000 ;
  171. FINSI ;
  172. SI ( Pe 'EGA' 100 ) ;
  173. LIUT = PROG 2.000 1.940 1.836 1.627 1.288 0.869 0.480 0.209
  174. 0.070 0.017 0.000 ;
  175. FINSI ;
  176. SI ( Pe 'EGA' 500 ) ;
  177. LIUT = PROG 2.000 2.000 1.998 1.965 1.702 0.947 0.242 0.023
  178. 0.001 0.000 0.000 ;
  179. FINSI ;
  180. SI ( Pe 'EGA' 1000 ) ;
  181. LIUT = PROG 2.000 2.000 2.000 1.985 1.841 0.951 0.154 0.001
  182. 0.000 0.000 0.000 ;
  183. FINSI ;
  184. SI ( Pe 'EGA' 1000000 ) ;
  185. LIUT = PROG 2.000 2.000 2.000 1.999 1.964 1.000 0.036 0.001
  186. 0.000 0.000 0.000 ;
  187. FINSI ;
  188.  
  189. EVOL4 = EVOL 'MANU' 'X' LIXT 'Uref(X,0)' LIUT ;
  190.  
  191. LIUC = IPOL LIXT LIX LIU ;
  192. NP = DIME LIXT ;
  193. ERR0 = 0. ;
  194. REPETER BLOC1 NP ;
  195. UCAL = EXTRAIRE LIUC &BLOC1 ;
  196. UREF = EXTRAIRE LIUT &BLOC1 ;
  197. ERR0 = ERR0 + ((UCAL-UREF)*(UCAL-UREF)) ;
  198. FIN BLOC1 ;
  199. ERR0 = ERR0/NP ;
  200. ERR0 = ERR0 '**' 0.5 ;
  201. *
  202. *
  203. *- POST-TRAITEMENT
  204. *
  205. *
  206. 'SI' GRAPH ;
  207. *
  208. * Maillage
  209. 'TRAC' DOMTOT 'TITR' 'Maillage' ;
  210. *
  211. * Vitesse
  212. UNCH = 'VECT' CHVIT 0.1 'UX' 'UY' 'ROUGE' ;
  213. 'TRAC' UNCH DOMTOT CNT1 'TITR' 'Vitesse transportante' ;
  214. *
  215. * Concentration
  216. 'TRAC' DOMTOT (RV1 . 'INCO' . 'CN') CNT1 'TITR' 'Concentration' ;
  217. *
  218. * Comparaison avec la solution analytique
  219. TAB1 = TABLE ;
  220. TAB1 . TITRE = TABLE ;
  221. TAB1. TITRE . 1 = 'MOT' 'CAST3M' ;
  222. TAB1. TITRE . 2 = 'MOT' 'REFERENCE' ;
  223. TAB1. 2 = 'MARQ LOSA NOLI' ;
  224. 'TITR' 'Comparaison CAST3M/Référence' ;
  225. 'DESS' (EVOL3 ET EVOL4) 'LEGE' TAB1
  226. 'MIMA' 'XBOR' -1.0 1.0 'YBOR' -1.0 3.0 ;
  227. *
  228. FINSI ;
  229. *
  230. * Arret si problème
  231. *
  232. SI ( ERR0 > tolera ) ;
  233. 'MESS' 'Erreur de ' err0 ' > ' tolera ;
  234. ERREUR 5 ;
  235. FINSI ;
  236. 'MESS' 'Erreur de ' err0 ' < ' tolera ;
  237. FIN ;
  238.  
  239.  
  240.  
  241.  
  242.  
  243.  
  244.  
  245.  
  246.  

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