Télécharger cc3d1.dgibi

Retour à la liste

Numérotation des lignes :

  1. * fichier : cc3d1.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4. **
  5. ** --- 08 Novembre 1999 ---
  6. **
  7. ** TEST CAVITE CUBIQUE
  8. **
  9. ** teste KBBT NS en 3D + le Bi CG cc2d1.dgibi
  10. ** teste les elements a pression continue sur des tetrahedres
  11. ** et pyramides P1-P1 et Q1-Q1 algorithme de projection
  12. ** stabilises en prenant DT > DT limite
  13. ** Teste aussi la creation de volumes automatique
  14. ** Si on ne veut que des tetrahedres activer option elem tri3
  15.  
  16. GRAPH = 'N' ;
  17. err1=5.e-2;
  18. DISCR='LINE' ;
  19. kpress='MSOMMET';
  20. COMPLET = FAUX ;
  21.  
  22. SI(COMPLET) ;
  23. ds1=0.05 ; ds2=0.08 ;
  24. ITMAX=20 ;
  25. SINON ;
  26. ds1=0.1 ; ds2=0.1 ;
  27. ITMAX=6;
  28. FINSI ;
  29.  
  30. option dime 2 elem qua4 ;
  31. *option dime 2 elem tri3 ;
  32.  
  33. opti isov suli ;
  34. p1= 0 0 ; p12=0.5 0. ; p2= 1 0 ;
  35.  
  36. ab=p1 d dini ds1 dfin ds2 p12 d dini ds2 dfin ds1 p2 ;
  37. ab12= p12 d dini ds1 dfin ds2 (0.5 0.5) d dini ds2
  38. dfin ds1 (0.5 1.) ;
  39. mt= ab trans dini ds1 dfin ds2 (0 0.5) trans dini ds2
  40. dfin ds1 (0 0.5) ;
  41. bc=cote 2 mt ;
  42. cd=cote 3 mt ;
  43. da=cote 4 mt ;
  44. ct=ab et bc et cd et da ;
  45. elim ct 1.e-3 ;
  46. mt=ab bc cd da daller ;
  47.  
  48. mt=orie mt ;
  49.  
  50. opti dime 3 elem cub8;
  51.  
  52. prof=0.1;
  53. mth=mt plus (0 0 prof) ;
  54. ab= ab12 plus (0 0 prof) ;
  55. oeil = 10 10 100 ;
  56. nph=2 ;
  57. cav=mt volu nph mth ;
  58. *trace oeil cav;
  59.  
  60. f1=face 1 cav ;
  61. f2=face 2 cav ;
  62. f3=face 3 cav ;
  63.  
  64. psup=cd trans nph (0 0 prof) ;
  65. psup=(CHAN TRI3 psup) ;
  66. ff3= (CHAN TRI3 f3) ;
  67. ELIM (PSUP ET FF3) 1.e-4 ;
  68. env3= f1 et f2 et ff3 ;
  69.  
  70. cav = volu env3 ;
  71.  
  72. cav= chan cav quaf ;
  73. ab=chan ab quaf ;
  74.  
  75. elim (cav et psup et ab et ff3 ) 1.e-5 ;
  76.  
  77. $mt=mode cav 'NAVIER_STOKES' DISCR ;
  78. $AB=mode AB 'NAVIER_STOKES' DISCR ;
  79. DOMA $mt IMPR ;
  80.  
  81. DT=2. ;
  82. MU=1. ;
  83. RO= 400. ;
  84.  
  85. prep1=doma $mt kpress;
  86. bcp=elem prep1 POI1 (lect 1 ) ;
  87.  
  88. rv= eqex 'OMEGA' 1. 'NITER' 1 ITMA ITMAX
  89. 'OPTI' 'EF' 'IMPL' KPRESS 'SUPG'
  90. ZONE $mt OPER NS 1. 'UN' (MU/RO) INCO 'UN'
  91. 'OPTI' 'EF' 'BDF2'
  92. ZONE $mt OPER DFDT 1. 'UNM' 'UNMM' DT 'UN' MU INCO 'UN'
  93. ;
  94.  
  95. rv=eqex rv
  96. CLIM
  97. UN UIMP ff3 0. UN VIMP ff3 0. UN WIMP ff3 0.
  98. UN WIMP (f1 et f2) 0.
  99. UN UIMP PSUP 1.
  100. ;
  101.  
  102.  
  103. rv.inco= table inco ;
  104. rv.inco.tn = kcht $mt scal sommet 5. ;
  105. rv.inco.un = kcht $mt vect sommet (1.e-5 1.e-5 1.e-5) ;
  106. rv.inco.'UNM' = kcht $mt vect sommet (1.e-5 1.e-5 1.e-5) ;
  107. rv.inco.'UNMM' = kcht $mt vect sommet (1.e-5 1.e-5 1.e-5) ;
  108.  
  109. rv.inco.pres = kcht $mt scal KPRESS 1.e-5 ;
  110.  
  111. rv.'METHINV'.TYPINV=3 ;
  112. rv.'METHINV'.IMPINV=1 ;
  113. rv.'METHINV'.NITMAX=1000;
  114. rv.'METHINV'.PRECOND=3 ;
  115. rv.'METHINV'.RESID =1.e-10 ;
  116.  
  117.  
  118. RVP= EQEX
  119. 'OPTI' 'EF' KPRESS
  120. ZONE $mt OPER KBBT (-1.) INCO 'UN' 'PRES'
  121. CLIM PRES TIMP bcp 0.
  122. ;
  123.  
  124. rvp.'METHINV'.TYPINV=2 ;
  125. rvp.'METHINV'.IMPINV=1 ;
  126. rvp.'METHINV'.NITMAX=300;
  127. rvp.'METHINV'.PRECOND=3 ;
  128. rvp.'METHINV'.RESID =1.e-8 ;
  129. rvp.'METHINV' . 'FCPRECT'=100 ;
  130. rvp.'METHINV' . 'FCPRECI'=100 ;
  131.  
  132. RV.'PROJ'= RVP ;
  133.  
  134. exec rv ;
  135.  
  136. srti=doma $AB 'MAILLAGE' ;
  137. evolV = EVOL 'CHPO' (rv.'INCO'.'UN') UX (srti ) ;
  138. evx= extr evolV 'ORDO' ;
  139. list evx ;
  140. evy=extr evolV 'ABSC' ;
  141. evolV= evol 'MANU' 'Vitesse' evx 'Hauteur' evy ;
  142.  
  143. si ('EGA' graph 'O' );
  144. TAB1=TABLE;
  145. TAB1.'TITRE'=TABLE ;
  146. TAB1 . 1 ='MARQ REGU ' ;
  147.  
  148. TAB1.'TITRE' . 1 = mot 'Composante_UX ' ;
  149. DESS evolV 'TITX' 'R (m)' 'TITY' 'V (m/s)' LEGE TAB1 ;
  150.  
  151. c1=vect (rv.inco.'UN') 0.3 ux uy uz jaune ;
  152. *trace c1 mt ;
  153. trace c1 cav ;
  154. trace (rv.inco.'PRESSION') mt ;
  155. finsi ;
  156.  
  157. lrr=prog
  158. 3.76000E-38 2.84554E-02 4.47889E-02 6.01943E-02 7.05196E-02
  159. 6.69563E-02 4.10891E-02 4.04425E-03 5.94186E-02 0.26219
  160. 1.0000 ;
  161.  
  162. evx = ABS evx ;
  163. list evx ;
  164. ER=SOMM( abs (evx - lrr) ) *0.1;
  165. mess ' Ecart sur MSOMMET TET4 PYR5 ' er ;
  166. Si ( er > 5.e-3) ; erreur 5 ; finsi ;
  167. FIN ;
  168.  
  169.  
  170.  
  171.  
  172.  
  173.  
  174.  
  175.  
  176.  
  177.  
  178.  
  179.  
  180.  
  181.  
  182.  
  183.  
  184.  
  185.  
  186.  
  187.  
  188.  
  189.  

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