Télécharger cc2d1.dgibi

Retour à la liste

Numérotation des lignes :

  1. * fichier : cc2d1.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4. **
  5. ** --- 10 Novembre 1999 ---
  6. **
  7. ** TEST CAVITE CARREE A PAROI DEFILANTE RE=400
  8. **
  9. ** teste LAPL KONV (CENTREE) KBBT pression continue
  10. ** Algorithme de projection
  11. ** Elements P1-P1 Q1-Q1 2D stables si DT >= DT limite
  12. ** DT pas trop grand non plus sinon pas de convergence
  13. ** TRI3 QUA4
  14. ** Algorithme VPI1
  15.  
  16.  
  17.  
  18. OPTION ISOV 'SULI' ;
  19. OPTION TRACE X;
  20. GRAPH = 'N' ;
  21. COMPLET = FAUX;
  22. *COMPLET = VRAI;
  23.  
  24.  
  25. Rey = 400. ;
  26. ds1=0.01 ; ds2=0.1 ;
  27.  
  28. Si COMPLET;
  29. ITMAX = 20 ;
  30. err1=1.e-4 ;
  31. Sinon ;
  32. ITMAX = 10 ;
  33. err1=9.e-2 ;
  34. Finsi ;
  35.  
  36.  
  37. KPRESS = 'MSOMMET';
  38. DISCR = 'LINE' ;
  39. TYPK = 'QUA4' ;
  40. KSUPG = 'CENTREE';
  41.  
  42. option dime 2 elem TYPK ;
  43.  
  44. p1= 0 0 ; p12=0.5 0. ; p2= 1 0 ;
  45.  
  46. ab=p1 d dini ds1 dfin ds2 p12 d dini ds2 dfin ds1 p2 ;
  47.  
  48. mt0 = ab trans dini ds1 dfin ds2 (0 0.5)
  49. trans dini ds2 dfin ds1 (0 0.5) ;
  50.  
  51. mt1= ab trans dini ds1 dfin ds2 (0 0.5);
  52.  
  53. opti elem TRI3 ;
  54.  
  55. mt2= (cote 3 mt1) trans dini ds2 dfin ds1 (0 0.5) ;
  56. mt2= inve mt2 ;
  57. bc=cote 2 mt0; cd=cote 3 mt0; da=cote 4 mt0;
  58. AA=bc moin (0.5 0.) ;
  59. BB=ab plus (0. 0.5) ;
  60.  
  61. elim (mt0 et mt1 et mt2) 1.e-5 ;
  62.  
  63. mt = mt1 et mt2 ;
  64. elim (AA et BB et mt) 1.e-5 ;
  65. cnt= cont mt ;
  66. mt= chan mt QUAF ;
  67. $mt=mode mt 'NAVIER_STOKES' DISCR ;
  68.  
  69. MU =1. ;
  70. RO= Rey ;
  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 1 ) ;
  76.  
  77.  
  78. Mmt= doma $mt 'MAILLAGE' ;
  79. doma $mt 'IMPR' ;
  80.  
  81. CD1= chan CD POI1 ;
  82. CD = elem CD1 (lect 2 pas 1 ((nbel cd1) - 1) );
  83.  
  84. RV= eqex 'OMEGA' 1. 'NITER' 1 ITMA itmax 'FIDT' 1
  85. 'OPTI' 'EF' 'IMPL' KSUPG KPRESS 'DIV2'
  86. ZONE $mt 'OPER' 'LAPN' (MU/RO) 'INCO' 'UN'
  87. ZONE $mt 'OPER' 'KONV' 1. 'UNM' MU DT 'INCO' 'UN'
  88. 'OPTI' 'EF' 'BDF2'
  89. ZONE $mt OPER DFDT 1. 'UNM' 'UNMM' DT 'UN' MU INCO 'UN'
  90. ;
  91.  
  92.  
  93. RV= eqex RV CLIM
  94. UN UIMP CD 1. UN VIMP CD 0. UN UIMP DA 0. UN VIMP DA 0.
  95. UN UIMP AB 0. UN VIMP AB 0. UN UIMP BC 0. UN VIMP BC 0. ;
  96.  
  97. rv.'METHINV'.TYPINV=3 ;
  98. rv.'METHINV'.IMPINV=0 ;
  99. rv.'METHINV'.NITMAX=400;
  100. rv.'METHINV'.PRECOND=3 ;
  101. rv.'METHINV'.RESID =1.e-8 ;
  102. rv. 'METHINV' . 'FCPRECT'=1 ;
  103. rv. 'METHINV' . 'FCPRECI'=1 ;
  104.  
  105. RVP= EQEX
  106. 'OPTI' 'EF' KPRESS
  107. ZONE $mt OPER KBBT (-1.) INCO 'UN' 'PRES'
  108. CLIM PRES TIMP bcp 0.
  109. ;
  110.  
  111. rvp.'METHINV'.TYPINV=2 ;
  112. rvp.'METHINV'.IMPINV=0 ;
  113. rvp.'METHINV'.NITMAX=300;
  114. rvp.'METHINV'.PRECOND=3 ;
  115. rvp.'METHINV'.RESID =1.e-8 ;
  116. rvp.'METHINV' . 'FCPRECT'=100 ;
  117. rvp.'METHINV' . 'FCPRECI'=100 ;
  118.  
  119. RV.'PROJ'= RVP ;
  120. RV.'ITMA' = ITMAX ;
  121.  
  122. rv.inco= table inco ;
  123. rv.inco.'UN' = kcht $mt vect sommet (1.e-5 1.e-5) ;
  124. rv.inco.'UNM' = kcht $mt vect sommet (1.e-5 1.e-5) ;
  125. rv.inco.'UNMM' = kcht $mt vect sommet (1.e-5 1.e-5) ;
  126. rv.inco.'PRES'= kcht $mt scal kpress 0.;
  127.  
  128. rv.'TYPROJ'='VPI1';
  129.  
  130. exec rv ;
  131. exec rv ;
  132.  
  133. AA = chan AA quaf ;
  134. BB = chan BB quaf ;
  135. elim (mt et aa et bb ) 1.e-3 ;
  136. $AA=mode AA 'NAVIER_STOKES' DISCR ;
  137. $BB=mode BB 'NAVIER_STOKES' DISCR ;
  138. srti=doma $AA 'MAILLAGE' ;
  139. srth=doma $BB 'MAILLAGE' ;
  140. evolV = EVOL 'CHPO' (rv.'INCO'.'UN') UX (srti ) ;
  141. evolH = EVOL 'CHPO' (rv.'INCO'.'UN') UY (srth ) ;
  142. evx='EXTR' evolV 'ORDO' ;
  143.  
  144. evy='EXTR' evolV 'ABSC' ;
  145. evolV= evol 'MANU' 'Vitesse' evx 'Hauteur' evy ;
  146. rv.'EVOLV'=evolV ;
  147. rv.'EVOLH'=evolH ;
  148.  
  149. si ('EGA' graph 'O' );
  150. TAB1=TABLE;
  151. TAB1.'TITRE'=TABLE ;
  152. TAB1 . 1 ='MARQ REGU ' ;
  153. TAB1.'TITRE' . 1 = mot 'Composante_UX ' ;
  154. DESS evolV 'TITX' 'R (m)' 'TITY' 'V (m/s)' LEGE TAB1 ;
  155.  
  156. TAB1.'TITRE' . 1 = mot 'Composante_UY ' ;
  157. DESS evolH 'TITX' 'R (m)' 'TITY' 'V (m/s)' LEGE TAB1 ;
  158. c1=vect (rv.inco.'UN') 0.3 ux uy jaune ;
  159. trace c1 mt ;
  160.  
  161. mtp= doma $mt 'MMAIL' ;
  162. trace (rv.inco.pression) mtp (cont mtp) ;
  163.  
  164. finsi ;
  165.  
  166. lrr = prog
  167. .35691E-38 1.83857E-02 3.87990E-02 6.18035E-02 8.78830E-02
  168. .11787 .15314 .19540 .24501 .29606
  169. .32867 .31443 .22963 .12329 9.19502E-03
  170. 8.18992E-02 .16068 .22291 .27060 .30865
  171. .34925 .41007 .50336 .62601 .76143
  172. .89040 1.0000;
  173.  
  174. lr = ABS( 'EXTRAIRE' (rv.'EVOLV') 'ABSC') ;
  175. list lr ;
  176.  
  177. Dgsrti = doma $AA 'XXDIAGSI' ;
  178. li= 'EXTR' (evol chpo Dgsrti srti) 'ORDO' ;
  179.  
  180. ER=( SOMM( ((lr - lrr)*(lr - lrr))*li) )** 0.5 ;
  181. mess ' Ecart Q1-Q1 P1-P1 projection 2D ' er ;
  182. Si ( er > err1 ) ; erreur 5 ; finsi ;
  183. FIN ;
  184.  
  185.  
  186.  
  187.  
  188.  
  189.  
  190.  
  191.  
  192.  
  193.  
  194.  
  195.  
  196.  
  197.  

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