Télécharger vahldavis.dgibi

Retour à la liste

Numérotation des lignes :

  1. *****************************************************
  2. ************************************************************************
  3. ************************************************************************
  4. * fichier : vahldavis.dgibi *
  5. ** modifie le 15/06/2014 passage EQPR -> EQEX *
  6. *****************************************************
  7. COMPLET = FAUX ;
  8. SI ( COMPLET ) ;
  9. d1 = 0.02 ;
  10. d2 = 0.1;
  11. NITER = 5000 ;
  12. CFL = 1.0 ;
  13. SINON ;
  14. d1 = 0.05 ;
  15. d2 = 0.5 ;
  16. NITER = 600 ;
  17. CFL = 1.0 ;
  18. FINSI ;
  19. GRAPH = FAUX ;
  20. KPRESS='CENTRE';
  21. DISCR ='MACRO';
  22. BETA=1.;
  23.  
  24.  
  25. ***************************************************************
  26. * CAVITE CARREE - VAHL-DAVIS *
  27. * A. Chene/H. Paillere TTMF Aout 1997 *
  28. ***************************************************************
  29.  
  30. ******************************
  31. *ESTIMATION DE LA CONVERGENCE*
  32. ******************************
  33.  
  34. DEBPROC CALCUL;
  35. ARGU RVX*'TABLE';
  36.  
  37. RV = RVX.'EQEX';
  38.  
  39. DD = RV.PASDETPS.'NUPASDT';
  40. NN = DD/5;
  41. LO = (DD-(5*NN)) EGA 0;
  42.  
  43. SI (LO);
  44.  
  45. UN = RV.INCO.'UN';
  46. UNM1 = RV.INCO.'UNM1';
  47.  
  48. unx = kcht $MT scal sommet (exco 'UX' un);
  49. unm1x = kcht $MT scal sommet (exco 'UX' unm1);
  50. uny = kcht $MT scal sommet (exco 'UY' un);
  51. unm1y = kcht $MT scal sommet (exco 'UY' unm1);
  52.  
  53. ERRX = KOPS unx '-' unm1x;
  54. ERRY = KOPS uny '-' unm1y;
  55. ELIX = MAXI ERRX 'ABS';
  56. ELIY = MAXI ERRY 'ABS';
  57. ELIX = (LOG (ELIX + 1.0E-10))/(LOG 10.);
  58. ELIY = (LOG (ELIY + 1.0E-10))/(LOG 10.);
  59. MESSAGE 'ITER' RV.PASDETPS.'NUPASDT' 'ERREUR LINF' ELIX ELIY;
  60. RV.INCO.'UNM1' = KCHT $MT 'VECT' 'SOMMET' (RV.INCO.'UN');
  61. IT = PROG RV.PASDETPS.'NUPASDT';
  62. ER = PROG ELIY;
  63. RV.INCO.'IT' = (RV.INCO.'IT') ET IT;
  64. RV.INCO.'ER' = (RV.INCO.'ER') ET ER;
  65.  
  66. FINSI;
  67. as2 ama1 = 'KOPS' 'MATRIK' ;
  68. FINPROC as2 ama1 ;
  69.  
  70. ************
  71. * MAILLAGE *
  72. ************
  73.  
  74. OPTI ELEM QUA8;
  75.  
  76. p1 = 0. 0.;
  77. p15 = 0.5 0.;
  78. p2 = 1. 0.;
  79. p25 = 1. 0.5;
  80. p3 = 1. 1.;
  81. p35 = 0.5 1.;
  82. p4 = 0. 1.;
  83. p45 = 0. 0.5;
  84.  
  85. bas = p1 d dini d1 dfin d2 p15 d dini d2 dfin d1 p2;
  86. cdro = p2 d dini d1 dfin d2 p25 d dini d2 dfin d1 p3;
  87. haut = p3 d dini d1 dfin d2 p35 d dini d2 dfin d1 p4;
  88. cgau = p4 d dini d1 dfin d2 p45 d dini d2 dfin d1 p1;
  89.  
  90. cnt = bas et cdro et haut et cgau;
  91. mt = bas cdro haut cgau daller;
  92.  
  93. ********
  94. * MODE *
  95. ********
  96.  
  97. Mmt= chan mt quaf ;
  98. $mt = MODE Mmt 'NAVIER_STOKES' DISCR ;
  99. doma $mt 'IMPR' ;
  100. mt = doma $mt 'MAILLAGE' ;
  101.  
  102. **************
  103. * PARAMETRES *
  104. **************
  105.  
  106. Pr = 0.71;
  107. Ra = 1.e6;
  108. Gr = Ra/Pr;
  109. NU = 1/(Gr**0.5);
  110. ALF= NU/Pr;
  111. gb = 0 -1;
  112. uref = 1;
  113.  
  114. ***********************
  115. * CREATION DES TABLES *
  116. ***********************
  117.  
  118. RV = EQEX $MT 'ITMA' NITER 'ALFA' CFL
  119. 'ZONE' $MT 'OPER' CALCUL
  120. 'OPTI' 'SUPG'
  121. 'ZONE' $MT 'OPER' 'NS' NU GB 'TN' 0.5 'INCO' 'UN'
  122. 'OPTI' 'SUPG' 'MMPG'
  123. 'ZONE' $MT 'OPER' 'TSCAL' ALF 'UN' 0. 'INCO' 'TN'
  124. 'OPTI' 'CENTREE'
  125. ZONE $MT 'OPER' 'DFDT' 1. 'UN' 'DELTAT' INCO 'UN'
  126. ZONE $MT 'OPER' 'DFDT' 1. 'TN' 'DELTAT' INCO 'TN'
  127. ;
  128.  
  129. RV = EQEX RV
  130. 'CLIM' 'UN' 'UIMP' cnt 0.
  131. 'UN' 'VIMP' cnt 0.
  132. 'TN' 'TIMP' cgau 1.
  133. 'TN' 'TIMP' cdro 0.;
  134.  
  135. RVP = EQEX 'OPTI' 'EF' KPRESS
  136. 'ZONE' $MT OPER KBBT -1. beta INCO 'UN' 'PRES'
  137. ;
  138.  
  139. rvp.'METHINV'.TYPINV=1 ;
  140. rvp.'METHINV'.IMPINV=0 ;
  141. rvp.'METHINV'.NITMAX=300;
  142. rvp.'METHINV'.PRECOND=3 ;
  143. rvp.'METHINV'.RESID =1.e-8 ;
  144. rvp.'METHINV' . 'FCPRECT'=100 ;
  145. rvp.'METHINV' . 'FCPRECI'=100 ;
  146.  
  147. RV.'PROJ' =RVP ;
  148.  
  149. RV.INCO = TABLE INCO;
  150. RV.INCO.'UN' = kcht $MT VECT SOMMET (0. 0.);
  151. RV.INCO.'PRES'= kcht $MT SCAL KPRESS 0.;
  152. corx = coor 1 mt;
  153. RV.INCO.'TN' = kcht $MT SCAL SOMMET (1.-corx);
  154.  
  155. RV.INCO.'UNM1' = kcht $MT VECT SOMMET (1.E-5 1.E-5);
  156. RV.INCO.'IT' = PROG 1;
  157. RV.INCO.'ER' = PROG 0.;
  158.  
  159. EXEC RV;
  160.  
  161. *************
  162. * RESULTATS *
  163. *************
  164.  
  165. Mcdro=chan cdro quaf ;
  166. elim (mt et cdro) 1.e-3 ;
  167. $cdro = mode Mcdro 'NAVIER_STOKES' DISCR;
  168. DY = DOMA $cdro 'VOLUME';
  169. DYT = SOMT DY;
  170. elim (Mmt et Mcdro)1.e-3;
  171.  
  172. GRADT = KOPS RV.INCO.'TN' 'GRAD' $MT;
  173. DTDX = KCHT $MT 'SCAL' 'CENTRE' (EXCO 'UX' GRADT);
  174. DTDX = ELNO $MT DTDX;
  175.  
  176. DTDXp = KCHT $CDRO 'SCAL' 'SOMMET' DTDX;
  177. DTDXp = KOPS DTDXp '*' (-1.);
  178.  
  179. EVOL1 = EVOL 'CHPO' DTDXp SCAL CDRO;
  180. EVOL1 = EVOL1 COUL ROUG ;
  181.  
  182. DTDXpc = NOEL $CDRO DTDXp;
  183. num = KOPS DTDXpc '*' DY;
  184. num = SOMT NUM;
  185. num = (-1)*NUM/DYT;
  186.  
  187. MESSAGE 'NUSSELT MOYEN' NUM;
  188. MESSAGE 'NUSSELT MAX' (MAXI DTDXpc);
  189. MESSAGE 'NUSSELT MIN' (MINI DTDXpc);
  190.  
  191. SI ( (MAXI DTDXp) < 12.00 ) ;
  192. ERREUR 5 ;
  193. FINSI ;
  194. SI ( ABS (NUM) < 7.55 ) ;
  195. ERREUR 5 ;
  196. FINSI ;
  197. SI ( (MINI RV.INCO.'ER') > -4.95 ) ;
  198. ERREUR 5 ;
  199. FINSI ;
  200.  
  201. SI ( GRAPH ) ;
  202. OPTI ISOV SULI;
  203. trace mt 'TITR' 'MAILLAGE' ;
  204. trac RV.INCO.'TN' mt (cont mt) 14 'TITR' 'TEMPERATURE' ;
  205.  
  206. unch = vect RV.INCO.'UN' 1. UX UY JAUNE;
  207. trace unch mt (cont mt) 'TITR' 'CHAMP DE VITESSE' ;
  208.  
  209. EVOL4 = EVOL 'MANU' 'ITERATIONS' (RV.INCO.'IT') 'LOG|E|inf'
  210. (RV.INCO.'ER') ;
  211. dess EVOL4 'XBOR' 0. 10000. 'YBOR' -10.0 0.0
  212. 'TITR' 'CONVERGENCE VERS LE STATIONNAIRE' ;
  213.  
  214. TAB1 = TABLE ;
  215. TAB1.1 = 'MARQ LOSA';
  216. DESS EVOL1 'XBOR' 0. 1. 'GRIL' TAB1 'TITR' 'NUSSELT A LA PAROI' ;
  217.  
  218. un = RV.INCO.'UN';
  219. un2=KOPS un '*' (Pr*(Gr**0.5));
  220. sw = (-1.) * (KOPS un 'ROT' $mt) ;
  221. rk = EQEX $mt 'OPTI' 'EF' 'IMPL'
  222. ZONE $mt OPER LAPN 1. INCO 'PSI'
  223. ZONE $mt OPER FIMP sw INCO 'PSI'
  224. 'CLIM' 'PSI' 'TIMP' (cgau et haut et bas et cdro) 0.;
  225. rk.'INCO'=table 'INCO' ;
  226. rk.'INCO'.'PSI'=kcht $mt scal sommet 0. ;
  227. EXEC rk ;
  228. psi=rk.'INCO'.'PSI';
  229. psi2=kops psi '*' (Pr*(Gr**0.5)) ;
  230. LISTE (MAXI PSI) ;
  231. LISTE (MINI PSI) ;
  232. trac psi mt CNT 14 'TITR' 'FONCTION DE COURANT' ;
  233. FINSI ;
  234.  
  235. FIN ;
  236.  
  237.  
  238.  
  239.  
  240.  
  241.  
  242.  
  243.  
  244.  
  245.  
  246.  
  247.  
  248.  

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