Télécharger cc2d3.dgibi

Retour à la liste

Numérotation des lignes :

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

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