Télécharger la plaquette Cast3M


fr en

Télécharger recirc.dgibi

Retour à la liste

Numérotation des lignes : 

  1. * fichier :  recirc.dgibi
  2. *
  3. *     Exemple test : recirculation dans une cavite semi circulaire
  4. *                    Eclt laminaire incompressible
  5. *                    On impose une vitesse tg sur la surface libre
  6. *                    V normale nulle sur le fond
  7. *
  8. *     Post  : debit sur le contour , fonction de courant , pression
  9. *     Teste vitesse normale imposee en 2D plan algorithme Taylor Hood
  10. *
  11. *     Voir aussi RECIRCP.DGIBI pour algorithme de projection
  12.  
  13. opti dime 2 elem QUA8 ;
  14. DISCR = 'QUAF';
  15. KPRES= 'CENTREP1' ;
  16.  
  17. GRAPH=FAUX ;
  18. *GRAPH=VRAI ;
  19. nbit = 5 ;
  20.  
  21. p0=0. 0. ;
  22. p1= -1. 0.;
  23. p2= 1. 0. ;
  24. p03= 0. -1. ;
  25. p3= 1. -1. ;
  26. p4= -1. -1. ;
  27.  
  28. n1=  5 ; n2= 5 ;
  29. ab= p1 d n1 p2 ;
  30. slib=  (chan ab POI1) ;
  31. nab= nbel slib ;
  32. slib= elem slib (lect 2 pas 1 (nab - 1));
  33.  
  34. bc= p2 c p0 n2 p03 c  p0 n2 p1;
  35. *bc= p2 d n2 p3 d n2 p4 d n2 p1 ;
  36.  
  37. cnt=ab et bc ;
  38.  
  39. mt= cnt surf ;
  40. Mmt=chan mt quaf ;
  41. Mbc=chan bc quaf ;
  42. elim (Mmt et Mbc) 1.e-5 ;
  43.  
  44. $mt = mode Mmt 'NAVIER_STOKES' DISCR;
  45. $bc = mode Mbc 'NAVIER_STOKES' DISCR;
  46.  
  47. * La cavité est fermée il faut imposer la pression en 1 point !
  48. prep1=doma $mt KPRES ;
  49. bcp=elem  prep1 POI1 (lect 1) ;
  50. *
  51.  
  52. *NU= 1. / 61. ;
  53. NU= 1.e-2 ;
  54. U0=1. ;
  55.  
  56. rv = eqex 'ITMA' 0 'NITER' nbit  'OMEGA' 0.8
  57.   OPTI EF 'CENTREE' 'IMPL' KPRES
  58.   'ZONE' $MT 'OPER' 'NS' 1. 'UN' NU  'INCO' 'UN'
  59.   'ZONE' $MT 'OPER' 'KBBT' 1.  'INCO' 'UN' 'PRES'
  60.   'ZONE' $bc 'OPER' 'VNIMP' $mt 'VN' 'INCO' 'UN' 'PRES'
  61.   'CLIM'
  62.          'UN' UIMP slib 1. 'UN' VIMP ab 0.
  63.          'PRES' TIMP bcp 0. ;
  64.  
  65. rv.'INCO'= TABLE 'INCO' ;
  66. rv.'INCO'.'UN'=kcht $mt vect sommet (0. 0.) ;
  67. rv.'INCO'.'PRES'=kcht $mt scal KPRES  0. ;
  68. rv.'INCO'.'VN'=kcht $bc scal sommet 0. ;
  69.  
  70. exec rv ;
  71.  
  72. un = rv.inco.'UN' ;
  73. ung= vect un 0.11  ux uy jaune ;
  74.  
  75. Si GRAPH ;trace ung Mmt ; Finsi ;
  76.  
  77. cnt1=chan cnt seg2 ;                                                    
  78. Mcnt1= chan cnt1 quaf;                                                  
  79. $cnt=model Mcnt1 'NAVIER_STOKES' line ;                                 
  80. q=dbit (rv.inco.'UN') $cnt ;                                            
  81. Mess ' Debit sur le contour ' q ;
  82. si ( q  > 1.5e-15 ) ; erreur 5 ; finsi ;
  83.  
  84. *
  85. * Calcul de la fonction de courant
  86. *
  87.  
  88. sw = 'KOPS' un 'ROT' $mt ;
  89. rk = 'EQEX' 'OPTI' 'EF' 'IMPL'
  90.      'ZONE' $mt 'OPER' 'LAPN' 1.0D0   'INCO' 'PSI'
  91.      'ZONE' $mt 'OPER' 'FIMP' sw      'INCO' 'PSI'
  92.      'CLIM' 'PSI' 'TIMP' cnt 0.0D0 ;
  93. rk . 'INCO' = 'TABLE' 'INCO' ;
  94. rk . 'INCO' . 'PSI' = 'KCHT' $mt 'SCAL' 'SOMMET' 0.0D0 ;
  95.  
  96. opti isov SULI ;
  97. exec rk ;
  98. psi=rk.inco.'PSI' ;
  99. Si GRAPH ;trace psi Mmt cnt ; Finsi ;
  100. mpsi=(mini psi);
  101. dpsi= abs (mpsi + .26820) ;
  102. mess ' MINI PSI = ' mpsi 'dpsi=' dpsi ;
  103. si ( dpsi > 1.e-5 ) ; erreur 5 ; finsi ;
  104.  
  105. *
  106. * Calcul de la pression
  107. *
  108.   pn= elno (rv.inco.'PRES') $mt KPRES;
  109. Si GRAPH ;  trace pn Mmt ; Finsi ;
  110.  
  111.  
  112. *
  113. *  Deuxieme cas : Meme probleme geometrie rectangulaire
  114. *                 MACRO CENTRE
  115. *
  116. *
  117. *
  118. *
  119. *
  120. *
  121.  
  122. DISCR = 'MACRO';
  123. KPRES= 'CENTRE' ;
  124.  
  125. p0=0. 0. ;
  126. p1= -1. 0.;
  127. p2= 1. 0. ;
  128. p03= 0. -1. ;
  129. p3= 1. -1. ;
  130. p4= -1. -1. ;
  131.  
  132. n1=  5 ; n2= 5 ;
  133. ab= p1 d n1 p2 ;
  134. slib=  (chan ab POI1) ;
  135. nab= nbel slib ;
  136. slib= elem slib (lect 2 pas 1 (nab - 1));
  137.  
  138. bc= p2 d n2 p3 d n2 p4 d n2 p1 ;
  139.  
  140. cnt=ab et bc ;
  141.  
  142. mt= cnt surf ;
  143. Mmt=chan mt quaf ;
  144. Mbc=chan bc quaf ;
  145. elim (Mmt et Mbc) 1.e-5 ;
  146.  
  147. $mt = mode Mmt 'NAVIER_STOKES' DISCR;
  148. $bc = mode Mbc 'NAVIER_STOKES' DISCR;
  149.  
  150. * La cavité est fermée il faut imposer la pression en 1 point !
  151. prep1=doma $mt KPRES ;
  152. bcp=elem  prep1 POI1 (lect 1) ;
  153. *
  154.  
  155. NU= 1.e-2 ;
  156. U0=1. ;
  157.  
  158. rv = eqex 'ITMA' 0 'NITER' nbit  'OMEGA' 0.8
  159.   OPTI EF 'CENTREE' 'IMPL' KPRES
  160.   'ZONE' $MT 'OPER' 'NS' 1. 'UN'  NU  'INCO' 'UN'
  161.   'ZONE' $MT 'OPER' 'KBBT' 1.  'INCO' 'UN' 'PRES'
  162.   'ZONE' $bc 'OPER' 'VNIMP' $mt 'VN' 'INCO' 'UN' 'PRES'
  163.   'CLIM'
  164.          'UN' UIMP slib 1. 'UN' VIMP ab 0.
  165.          'PRES' TIMP bcp 0. ;
  166.  
  167. rv.'INCO'= TABLE 'INCO' ;
  168. rv.'INCO'.'UN'=kcht $mt vect sommet (0. 0.) ;
  169. rv.'INCO'.'PRES'=kcht $mt scal KPRES 0. ;
  170. rv.'INCO'.'VN'=0. ;
  171.  
  172. exec rv ;
  173.  
  174.  
  175. un = rv.inco.'UN' ;
  176. ung= vect un 0.11  ux uy jaune ;
  177.  
  178. Si GRAPH ;trace ung Mmt ; Finsi ;
  179.  
  180. cnt1=chan cnt seg2 ;                                                    
  181. Mcnt1= chan cnt1 quaf;                                                  
  182. $cnt=model Mcnt1 'NAVIER_STOKES' line ;                                 
  183. q=dbit (rv.inco.'UN') $cnt ;                                            
  184. Mess ' Debit sur le contour ' q ;
  185. si ( q  > 1.5e-15 ) ; erreur 5 ; finsi ;
  186.  
  187. sw = 'KOPS' un 'ROT' $mt ;
  188. rk = 'EQEX' 'OPTI' 'EF' 'IMPL'
  189.      'ZONE' $mt 'OPER' 'LAPN' 1.0D0   'INCO' 'PSI'
  190.      'ZONE' $mt 'OPER' 'FIMP' sw      'INCO' 'PSI'
  191.      'CLIM' 'PSI' 'TIMP' cnt 0.0D0 ;
  192. rk . 'INCO' = 'TABLE' 'INCO' ;
  193. rk . 'INCO' . 'PSI' = 'KCHT' $mt 'SCAL' 'SOMMET' 0.0D0 ;
  194.  
  195. opti isov SULI ;
  196. exec rk ;
  197. psi=rk.inco.'PSI' ;
  198. mess ' MINI PSI = ' (mini psi) ;
  199. Si GRAPH ;trace psi Mmt cnt ; Finsi ;
  200. mpsi=(mini psi);
  201. dpsi= abs (mpsi + .22) ;
  202. mess ' MINI PSI = ' mpsi 'dpsi=' dpsi ;
  203. si ( dpsi > 1.e-3 ) ; erreur 5 ; finsi ;
  204.  
  205.  
  206. *
  207. * Calcul de la pression
  208. *
  209.   pn= elno (rv.inco.'PRES') $mt KPRES;
  210. Si GRAPH ;  trace pn Mmt ; Finsi ;
  211. FIN ;
  212.  
  213.  
  214.  
  215.  
  216.  
  217.  
  218.  
  219.  
  220.  
  221.  
  222.  

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