Télécharger ccar5.dgibi

Retour à la liste

Numérotation des lignes :

  1. * fichier : ccar5.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4. **
  5. ** --- 6 JUIN 1999 ---
  6. **
  7. ** TEST CAVITE CARRE A PAROI DEFILANTE RE=400
  8. **
  9. ** Comparaison des algorithmes Implicite, projection
  10. ** et semi explicite.
  11. ** teste LAPL KONV NS (CENTREE) KBBT
  12. ** Formulation MACRO CENTRE (Stabilisation) ELEMENT QUA9
  13. ** La comparaison fine se fait avec COMPLET=VRAI
  14. ** prevoir espace et temps cpu
  15.  
  16. OPTION ISOV 'SULI' ;
  17. *OPTION TRACE PS;
  18. GRAPH = 'N' ;
  19. GRAPHI= 'N' ;
  20. COMPLET = FAUX ;
  21. *COMPLET = VRAI ;
  22. Rey = 400.;
  23.  
  24. SI ( COMPLET ) ;
  25.  
  26. ds1=0.02; ds2=0.1 ;
  27. ds1=0.01; ds2=0.08;
  28. ITMAX = 40 ;
  29. ITMAI = 30 ;
  30. ITMAE = 5000 ;
  31. KSUPG = 'CENTREE';
  32. SINON ;
  33.  
  34. err1=1.e-4 ;
  35. ds1=0.2; ds2=0.2 ;
  36. ITMAX = 30 ;
  37. ITMAI = 30 ;
  38. ITMAE = 2000;
  39. ITMAX = 10;
  40. ITMAI = 10;
  41. ITMAE = 100;
  42. KSUPG = 'SUPG';
  43.  
  44. FINSI ;
  45.  
  46. KPRESS = 'CENTRE' ;
  47. DISCR = 'LINE' ;
  48. TYPK = 'QUA8' ;
  49. ALGO = 'SPLT' ;
  50.  
  51. DEBPROC TEST KPRESS*MOT TYPK*MOT DISCR*MOT GRAPH*MOT ALGO*MOT
  52. ITMAX*ENTIER ITMAI*ENTIER ITMAE*ENTIER ;
  53. option dime 2 elem TYPK ;
  54.  
  55. p1= 0 0 ; p12=0.5 0. ; p2= 1 0 ;
  56.  
  57. ab=p1 d dini ds1 dfin ds2 p12 d dini ds2 dfin ds1 p2 ;
  58. mt= ab trans dini ds1 dfin ds2 (0 0.5) trans dini ds2
  59. dfin ds1 (0 0.5) ;
  60. bc=cote 2 mt ; cd=cote 3 mt ; da=cote 4 mt ;
  61. AA=bc moin (0.5 0.) ;
  62. BB=ab plus (0. 0.5) ;
  63. elim (AA et BB et mt) 1.e-3 ;
  64.  
  65. mt= chan mt QUAF ;
  66. $mt=mode mt 'NAVIER_STOKES' DISCR ;
  67.  
  68. MU =1. ;
  69. RO= Rey ;
  70. EPSS=-1.e-8 ;
  71. DT=1. ;
  72.  
  73. * La cavité est fermée il faut imposer la pression en 1 point !
  74. prep1=doma $mt kpress;
  75. bcp=elem prep1 POI1 (lect 15);
  76. bcp = bcp coul rouge ;
  77. Mmt= doma $mt 'MAILLAGE' ;
  78. doma $mt 'IMPR' ;
  79. Si (EGA GRAPH 'O') ; trace (Mmt et bcp); Finsi ;
  80. *
  81.  
  82. Mab = chan quaf ab ;
  83. Mcd = chan quaf cd ;
  84. Mda = chan quaf da ;
  85. Mbc = chan quaf bc ;
  86. elim (Mmt et Mab et Mcd et Mda et Mbc) 1.e-5 ;
  87. $ab = model Mab 'NAVIER_STOKES' DISCR ;
  88. $cd = model Mcd 'NAVIER_STOKES' DISCR ;
  89. $da = model Mda 'NAVIER_STOKES' DISCR ;
  90. $bc = model Mbc 'NAVIER_STOKES' DISCR ;
  91. ab = doma $ab maillage ;
  92. cd = doma $cd maillage ;
  93. da = doma $da maillage ;
  94. bc = doma $bc maillage ;
  95.  
  96. Si (EGA ALGO 'EXPL');
  97.  
  98. RV= eqex 'OMEGA' 1. 'NITER' 1 ITMA itmax 'FIDT' 1
  99. 'OPTI' 'EFM1' KSUPG
  100. ZONE $mt OPER NS (MU/RO) INCO 'UN'
  101. 'OPTI' 'EFM1' 'CENTREE'
  102. ZONE $mt OPER DFDT 1. 'UN' 'DELTAT' INCO 'UN' ;
  103.  
  104. Sinon ;
  105.  
  106. RV= eqex 'OMEGA' 1. 'NITER' 1 ITMA itmax 'FIDT' 1
  107. 'OPTI' 'EF' 'IMPL' KSUPG
  108. * ZONE $mt OPER NS (MU/RO) INCO 'UN'
  109. ZONE $mt 'OPER' 'LAPN' MU 'INCO' 'UN'
  110. ZONE $mt 'OPER' 'KONV' RO 'UN' MU DT 'INCO' 'UN'
  111. 'OPTI' 'EFM1' 'CENTREE'
  112. ZONE $mt OPER DFDT RO 'UN' DT INCO 'UN'
  113. ;
  114.  
  115. Finsi ;
  116.  
  117. RV= eqex RV CLIM
  118. UN UIMP CD 1. UN VIMP CD 0. UN UIMP DA 0. UN VIMP DA 0.
  119. UN UIMP AB 0. UN VIMP AB 0. UN UIMP BC 0. UN VIMP BC 0. ;
  120.  
  121. rv.'METHINV'.TYPINV=3 ;
  122. rv.'METHINV'.IMPINV=0 ;
  123. rv.'METHINV'.NITMAX=400;
  124. rv.'METHINV'.PRECOND=3 ;
  125. rv.'METHINV'.RESID =1.e-8 ;
  126. rv. 'METHINV' . 'FCPRECT'=1 ;
  127. rv. 'METHINV' . 'FCPRECI'=1 ;
  128.  
  129. 'SI' (EGA ALGO 'IMPL') ;
  130.  
  131. RV.'ITMA' = ITMAI ;
  132. RV.'OMEGA' = 1. ;
  133.  
  134. RV= eqex RV
  135. 'OPTI' 'EF' 'IMPL' KSUPG KPRESS
  136. ZONE $mt OPER KBBT 1. INCO 'UN' 'PRES'
  137. CLIM PRES TIMP bcp 0. ;
  138.  
  139. 'SINON' ;
  140.  
  141. Si (EGA ALGO 'EXPL') ; betastab = 1. ;
  142. Sinon ; betastab= 1.e2 ;
  143. Finsi ;
  144. Si (EGA TYPK 'TRI6') ; betastab = 1.e5 ; Finsi ;
  145.  
  146. RVP= EQEX
  147. 'OPTI' 'EFM1' KPRESS
  148. ZONE $mt OPER KBBT -1. betastab INCO 'UN' 'PRES'
  149. CLIM PRES TIMP bcp 0. ;
  150.  
  151. rvp.'METHINV'.TYPINV=2 ;
  152. rvp.'METHINV'.TYPINV=1 ;
  153. rvp.'METHINV'.IMPINV=0 ;
  154. rvp.'METHINV'.NITMAX=300;
  155. rvp.'METHINV'.PRECOND=3 ;
  156. rvp.'METHINV'.RESID =1.e-8 ;
  157. rvp.'METHINV' . 'FCPRECT'=100 ;
  158. rvp.'METHINV' . 'FCPRECI'=100 ;
  159.  
  160. si (EGA ALGO 'SPLT') ;
  161. RV.'PROJ'= RVP ;
  162. RV.'ITMA' = ITMAX ;
  163. Finsi ;
  164.  
  165. si (EGA ALGO 'EXPL') ;
  166. RV.'PROJ' = RVP ;
  167. RV.'TYPROJ'='PSCT';
  168. RV.'ITMA' = ITMAE ;
  169. RV.'ALFA' = 0.5 ;
  170. RV.'FIDT' = 50 ;
  171. Finsi ;
  172.  
  173. 'FINSI' ;
  174.  
  175. rv.inco= table inco ;
  176. rv.inco.un = kcht $mt vect sommet (1.e-5 1.e-5) ;
  177. rv.inco.'PRES'= kcht $mt scal kpress 0.;
  178.  
  179. exec rv ;
  180. exec rv ;
  181.  
  182. AA = chan AA quaf ;
  183. BB = chan BB quaf ;
  184. elim (mt et aa et bb ) 1.e-3 ;
  185. $AA=mode AA 'NAVIER_STOKES' DISCR ;
  186. $BB=mode BB 'NAVIER_STOKES' DISCR ;
  187. srti=doma $AA 'MAILLAGE' ;
  188. srth=doma $BB 'MAILLAGE' ;
  189. evolV = EVOL 'CHPO' (rv.'INCO'.'UN') UX (srti ) ;
  190. evolH = EVOL 'CHPO' (rv.'INCO'.'UN') UY (srth ) ;
  191. evx='EXTR' evolV 'ORDO' ;
  192. list evx ;
  193. evy='EXTR' evolV 'ABSC' ;
  194. evolV= evol 'MANU' 'Vitesse' evx 'Hauteur' evy ;
  195. rv.'EVOLV'=evolV ;
  196. rv.'EVOLH'=evolH ;
  197.  
  198. si ('EGA' graph 'O' );
  199. TAB1=TABLE;
  200. TAB1.'TITRE'=TABLE ;
  201. TAB1 . 1 ='MARQ REGU ' ;
  202. TAB1.'TITRE' . 1 = mot 'Composante_UX ' ;
  203. DESS evolV 'TITX' 'R (m)' 'TITY' 'V (m/s)' LEGE TAB1 ;
  204.  
  205. TAB1.'TITRE' . 1 = mot 'Composante_UY ' ;
  206. DESS evolH 'TITX' 'R (m)' 'TITY' 'V (m/s)' LEGE TAB1 ;
  207. c1=vect (rv.inco.'UN') 0.3 ux uy jaune ;
  208. trace c1 mt ;
  209.  
  210. Si (EGA ALGO 'IMPL') ;
  211. pn= elno $mt (kcht $mt scal centre ( nomc scal (rv.inco.'PRES')));
  212. sinon ;
  213. pn= elno $mt (kcht $mt scal centre ( nomc scal (rv.inco.pression)));
  214. finsi ;
  215. trace pn mt ;
  216.  
  217. mess ' MAX P ' (maxi pn) 'MIN P ' (mini pn) ;
  218. finsi ;
  219. FINPROC RV ;
  220.  
  221.  
  222.  
  223. RV1= test 'CENTRE' TYPK 'MACRO' GRAPH 'SPLT' ITMAX ITMAI ITMAE ;
  224.  
  225. evv= (rv1.'EVOLV');
  226. evh= (rv1.'EVOLH');
  227.  
  228. SI ( NON COMPLET ) ;
  229.  
  230. lr = 'EXTRAIRE' (rv1.'EVOLV') 'ABSC' ;
  231. list lr ;
  232.  
  233. lrr = prog
  234. -8.26684E-38 -2.74863E-02 -5.08159E-02 -7.49178E-02 -9.42076E-02
  235. -9.99882E-02 -8.10852E-02 -4.32674E-02 7.12863E-03 5.44891E-02
  236. 0.10073 0.30942 1.0000;
  237. ER=SOMM( abs (lr - lrr) ) ;
  238. mess ' Ecart sur CENTRE QUA8 MACRO Projection ' er ;
  239. Si ( er > err1 ) ; erreur 5 ; finsi ;
  240. FINSI ;
  241.  
  242.  
  243. RV2= test 'CENTRE' TYPK 'MACRO' GRAPH 'IMPL' ITMAX ITMAI ITMAE ;
  244.  
  245. evv= evv et (rv2.'EVOLV');
  246. evh= evh et (rv2.'EVOLH');
  247.  
  248. SI ( NON COMPLET ) ;
  249.  
  250. lr = 'EXTRAIRE' (rv2.'EVOLV') 'ABSC' ;
  251. list lr ;
  252.  
  253. lrr=prog
  254. 0.0000 -3.10947E-02 -5.38085E-02 -7.39196E-02 -8.54276E-02
  255. -8.30593E-02 -6.45006E-02 -3.23991E-02 3.13243E-03 3.90764E-02
  256. 8.72317E-02 0.29477 1.0000;
  257.  
  258. ER=SOMM( abs (lr - lrr) ) ;
  259. mess ' Ecart sur CENTRE QUA8 MACRO Implicite ' er ;
  260. Si ( er > err1 ) ; erreur 5 ; finsi ;
  261. FINSI ;
  262.  
  263.  
  264. RV3= test 'CENTRE' TYPK 'MACRO' GRAPH 'EXPL' ITMAX ITMAI ITMAE ;
  265. list RV3.'TYPROJ' ;
  266.  
  267. evv= evv et (rv3.'EVOLV');
  268. evh= evh et (rv3.'EVOLH');
  269.  
  270. SI ( NON COMPLET ) ;
  271.  
  272. lr = 'EXTRAIRE' (rv3.'EVOLV') 'ABSC' ;
  273. list lr ;
  274.  
  275. lrr=prog
  276. 3.29048E-37 -1.67015E-02 -2.71125E-02 -3.64044E-02 -4.66237E-02
  277. -5.61821E-02 -5.96575E-02 -5.28560E-02 -3.42689E-02 -2.85057E-03
  278. 5.39080E-02 0.27875 1.0000;
  279.  
  280. ER=SOMM( abs (lr - lrr) ) ;
  281. mess ' Ecart sur CENTRE TRI6 MACRO Projection ' er ;
  282. Si ( er > err1 ) ; erreur 5 ; finsi ;
  283. FINSI ;
  284. *FIN ;
  285.  
  286. si ('EGA' graphi 'O' );
  287.  
  288. * Solution de Ghia et al (1982) :
  289.  
  290. y = 'PROG' 0. 0.0547 0.0625 0.0703 0.1016 0.1719 0.2813 0.4531 0.5
  291. 0.6172 0.7344 0.8516 0.9531 0.9609 0.9688 0.9766 1.;
  292.  
  293. 'SI' (Rey 'EGA' 100.) ;
  294. u = 'PROG' 0. -0.03717 -0.04192 -0.04775 -0.06434 -0.10150 -0.15662
  295. -0.2109 -0.20581 -0.13641 0.00332 0.23151 0.68717 0.73722
  296. 0.78871 0.84123 1.;
  297. 'FINSI' ;
  298.  
  299. 'SI' (Rey 'EGA' 400.) ;
  300. u = 'PROG' 0. -0.08186 -0.09266 -0.10338 -0.14612 -0.24299 -0.32726
  301. -0.17119 -0.11477 0.02135 0.16256 0.29093 0.55892 0.61756
  302. 0.68439 0.75837 1.;
  303. 'FINSI' ;
  304.  
  305. 'SI' (Rey 'EGA' 1000.) ;
  306. u = 'PROG' 0. -0.18109 -0.20196 -0.22220 -0.29730 -0.38289 -0.27805
  307. -0.10648 -0.06080 0.05702 0.18719 0.33304 0.46604 0.51117
  308. 0.57492 0.65928 1.;
  309. 'FINSI' ;
  310.  
  311. evolug = 'EVOL' 'MANU' u y;
  312.  
  313. TAB1=TABLE;
  314. TAB1.'TITRE'=TABLE ;
  315. TAB1 . 1 ='MARQ CROI REGU ' ;
  316. TAB1 . 2 ='MARQ TRIA REGU ' ;
  317. TAB1 . 3 ='MARQ LOSA REGU ' ;
  318. TAB1 . 4 ='TIRR MARQ PLUS REGU';
  319. TAB1.'TITRE' . 1 = mot ' Proj ' ;
  320. TAB1.'TITRE' . 2 = mot ' Impl ' ;
  321. TAB1.'TITRE' . 3 = mot ' 1/2 Expl' ;
  322. TAB1.'TITRE' . 4 = MOT 'V Ghia et al (1982)';
  323.  
  324. evv = evv et evolug ;
  325.  
  326. DESS evv 'TITX' 'R (m)' 'TITY' 'V (m/s)' LEGE TAB1
  327. DATE MIMA TITR
  328. ' Cavite carree (Qua9 Macro stab comp ux) ; Comparaison Algorithmes ';
  329.  
  330. * Solution de Ghia et al (1982) :
  331.  
  332. x = 'PROG' 0. 0.0625 0.0703 0.0781 0.0938 0.1563 0.2266 0.2344 0.5
  333. 0.8047 0.8594 0.9063 0.9453 0.9531 0.9609 0.9688 1.;
  334.  
  335. 'SI' (Rey 'EGA' 100.);
  336. v = 'PROG' 0. 0.09233 0.10091 0.10890 0.12317 0.16077 0.17507 0.17527
  337. 0.05454 -0.24533 -0.22445 -0.16914 -0.10313 -0.08864
  338. -0.07391 -0.05906 0.;
  339. 'FINSI' ;
  340.  
  341. 'SI' (Rey 'EGA' 400.);
  342. v = 'PROG' 0. 0.18360 0.19713 0.20920 0.22965 0.28124 0.30203 0.30174
  343. 0.05186 -0.38598 -0.44993 -0.33827 -0.22847 -0.19254
  344. -0.15663 -0.12146 0.;
  345. 'FINSI' ;
  346. 'SI' (Rey 'EGA' 1000.);
  347. v = 'PROG' 0. 0.27485 0.29012 0.30353 0.32627 0.37095 0.33075 0.32235
  348. 0.02526 -0.31966 -0.42665 -0.51550 -0.39188 -0.33714
  349. -0.27669 -0.21388 0.;
  350. 'FINSI' ;
  351.  
  352. evolvg = 'EVOL' 'MANU' x v;
  353.  
  354. evh = evh et evolvg ;
  355.  
  356.  
  357. DESS evh 'TITX' 'R (m)' 'TITY' 'V (m/s)' LEGE TAB1
  358. DATE MIMA TITR
  359. ' Cavite carree (Qua9 Macro stab comp uy) ; Comparaison Algorithmes ';
  360. finsi ;
  361.  
  362. *Si complet ;
  363. *opti sauv 'ccar5.dgibi.sauv' ;
  364. *sauv ;
  365. *finsi ;
  366. FIN ;
  367.  
  368.  
  369.  
  370.  
  371.  
  372.  
  373.  
  374.  
  375.  
  376.  
  377.  
  378.  
  379.  
  380.  
  381.  

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