Télécharger coude.dgibi

Retour à la liste

Numérotation des lignes :

  1. * fichier : coude.dgibi
  2. opti dime 2 elem QUA8 ;
  3. *opti dime 2 elem TRI6 ;
  4. GRAPH=VRAI ;
  5. GRAPH=FAUX ;
  6.  
  7. DISCR='MACRO' ;
  8. *DISCR='QUAF';
  9. *DISCR='LINE';
  10. KPRES='CENTREP1' ;
  11. *KPRES='CENTRE' ;
  12.  
  13.  
  14. *R=0.5 ; nl = 4 ;n1= 5 ; lgd1=1.; lgd2=2.; n2=5 ;
  15. R=0.5 ; nl = 2 ;n1= 2 ; lgd1=1.; lgd2=2.; n2=2 ;
  16. P0=0 0 ; p1= R 0 ;
  17. p2 = p1 tour 90 p0 ;
  18. p3 = p2 tour 90 p0 ;
  19. p4 = p3 tour 90 p0 ;
  20. mct= p1 c p0 n1 p2 c p0 n1 p3 c p0 n1 p4 c p0 n1 p1;
  21. ct1= p1 c p0 n1 p2;
  22. ct2= p2 c p0 n1 p3;
  23. ct3= p3 c p0 n1 p4;
  24. ct4= p4 c p0 n1 p1;
  25.  
  26. *ct2= p2 d n1 p3;
  27. *ct4= p4 d n1 p1;
  28.  
  29. *mct= p1 c p0 n1 p2 d n1 p0 d n1 p1;
  30. *mc= surf mct ;
  31. mc= daller ct1 ct2 ct3 ct4 ;
  32. *trace mc ;
  33.  
  34. opti dime 3 elem seg2 ;
  35.  
  36. q1= p0 plus (0 0 lgd1);
  37. ax1 = q1 plus (1 0 0) ;
  38. ax2 = q1 plus (1 1 0) ;
  39. q2= q1 tour 90. ax1 ax2 ;
  40. q3= q2 plus (lgd2 0 0);
  41.  
  42. lg= p0 d nl q1 c ax1 n2 q2 d nl q3 ;
  43. *lg= p0 d nl q1 ;
  44.  
  45. opti elem CU20 ;
  46.  
  47. *trace lg ;
  48. lg1= chan lg 'POI1' ;
  49. nbp= (nbno lg1) - 2 ;
  50. pi1= p0 ;
  51. si1= mc ;
  52. pi2= poin 2 lg1 ;
  53. vh1=pi2 moins pi1 ;
  54. si2= si1 plus vh1 ;
  55. Entree = si1 ;
  56. vc=si1 volu 2 si2;
  57. *trace cache vc ;
  58.  
  59. repeter bloc1 nbp ;
  60. pi1=pi2;
  61. si1=si2 ;
  62. pi2= poin (&bloc1 + 2) lg1 ;
  63. vh2=pi2 moins pi1 ;
  64. vhr=vh2 ;
  65.  
  66. axe =vh1 pvec vhr ;
  67. alfa= atg (norm(axe / (vh1 psca vhr) ));
  68. Si ( alfa >EG 1.) ;
  69. si2= si1 tour alfa pi1 (pi1 plus axe) ;
  70. Sinon ;
  71. si2= si1 ;
  72. Finsi ;
  73. si2=si2 plus vh2 ;
  74. vh1=vh2 ;
  75. vc = vc et (si1 volu 2 si2);
  76. fin bloc1 ;
  77. Sortie = si2 ;
  78. oeil = 0. -1.e5 -100. ;
  79. *trace cache vc oeil ;
  80. Sortie= Sortie coul rouge ;
  81. Entree= Entree coul rouge ;
  82.  
  83. Mvc = chan vc QUAF ;
  84. Mentree = chan entree QUAF ;
  85. Msortie = chan sortie QUAF ;
  86.  
  87. ELIM (Mvc et Mentree et Msortie) 1.e-5 ;
  88.  
  89. $vc = mode Mvc 'NAVIER_STOKES' DISCR ;
  90. $vcm= mode Mvc 'NAVIER_STOKES' 'MACRO';
  91.  
  92. $sortie = mode Msortie 'NAVIER_STOKES' DISCR ;
  93.  
  94. $entree = mode Mentree 'NAVIER_STOKES' DISCR ;
  95. entree = doma $entree maillage ;
  96. Menvlp = doma $vc 'ENVELOPPE' ;
  97. Mparoi = diff Menvlp (Mentree et Msortie) ;
  98.  
  99. doma $vc 'IMPR' ;
  100.  
  101. * On extrait le champ des normales de la surface d'entree
  102. * en prenant garde a l'orientation des normales
  103. nentr= doma $entree 'NORMALEV';
  104.  
  105. Mcte = doma $entree 'ENVELOPPE';
  106. $cte= mode Mcte 'NAVIER_STOKES' DISCR ;
  107. cte=doma $cte maillage ;
  108.  
  109. * On extrait les composantes nx ny etc
  110. nx= exco 'UX' nentr;
  111. ny= exco 'UY' nentr;
  112. nz= exco 'UZ' nentr;
  113. * on cree les componsantes cartesiennes du champ de vitesse
  114. * qui seront imposees en condition limite type Dirichlet
  115. U0=1. ;
  116. usx=u0*nx ;
  117. usy=u0*ny ;
  118. usz=u0*nz ;
  119.  
  120. $envlp = mode Menvlp 'NAVIER_STOKES' DISCR ;
  121. envlp=doma $envlp maillage ;
  122. $paroi = mode Mparoi 'NAVIER_STOKES' DISCR ;
  123. $parom = mode Mparoi 'NAVIER_STOKES' 'MACRO';
  124. paroi=doma $paroi maillage ;
  125.  
  126. NU=1.e-5 ;
  127. NUT=5.e-3;
  128. YP=5.e-3 ;
  129. DT=1. ;
  130.  
  131. rv= eqex 'ITMA' 10 'NITER' 1 'OMEGA' 1.
  132. OPTI EF 'SUPG' 'IMPL' KPRES
  133. 'ZONE' $paroi 'OPER' 'FPU' 1. 'UN' NU 'UET' YP 'INCO' 'UN'
  134. 'ZONE' $vc 'OPER' 'NS' 1. 'UN' NUT 'INCO' 'UN'
  135. * 'ZONE' $vc 'OPER' 'KBBT' -1. 'INCO' 'UN' 'PRES'
  136. * 'ZONE' $paroi 'OPER' 'VNIMP' 0. 'INCO' 'UN'
  137. OPTI EF 'CENTREE'
  138. 'ZONE' $vc 'OPER' 'DFDT' 1. 'UN' DT 'INCO' 'UN'
  139. 'CLIM'
  140. 'UN' UIMP entree usx 'UN' VIMP entree usy 'UN' WIMP entree usz
  141. * 'UN' UIMP paroi 0. 'UN' VIMP paroi 0. 'UN' WIMP paroi 0.
  142. ;
  143.  
  144. rv.'METHINV'.TYPINV=3 ;
  145. rv.'METHINV'.IMPINV=0 ;
  146. rv.'METHINV'.NITMAX=100;
  147. rv.'METHINV'.PRECOND=3 ;
  148. rv.'METHINV'.RESID =1.e-10 ;
  149.  
  150. rv.'INCO'= TABLE 'INCO' ;
  151. rv.'INCO'.'UN'=kcht $vc vect sommet (0. 0. 0.) ;
  152. rv.'INCO'.'PRES'=kcht $vc scal KPRES 0. ;
  153.  
  154. RVP= EQEX 'OPTI' 'EF' KPRES
  155. 'ZONE' $vc 'OPER' 'KBBT' -1. 'INCO' 'UN' 'PRES'
  156. 'ZONE' $paroi 'OPER' 'VNIMP' $vc 0. 'INCO' 'UN' 'PRES'
  157. ;
  158.  
  159. rvp.'METHINV'.TYPINV=3 ;
  160. rvp.'METHINV'.IMPINV=0 ;
  161. rvp.'METHINV'.NITMAX=100;
  162. rvp.'METHINV'.PRECOND=3 ;
  163. rvp.'METHINV'.RESID =1.e-10 ;
  164.  
  165.  
  166. RV.'PROJ'= RVP ;
  167.  
  168. exec rv ;
  169.  
  170. un=rv.inco.'UN' ;
  171. un1=redu un paroi ;
  172. ung=vect un1 0.1 ux uy uz jaune ;
  173.  
  174. Si GRAPH ;
  175. trace ung paroi oeil ;
  176. Finsi;
  177.  
  178. q= dbit un $envlp ;
  179. list q ;
  180. qe=dbit un $entree ;
  181. qs=dbit un $sortie ;
  182. qp=dbit un $paroi ;
  183. mess 'qe=' qe 'qs=' qs 'qp=' qp 'qs+qp=' (qs+qp) ;
  184.  
  185.  
  186. yplus=yp*(rv.inco.'UET')/nu ;
  187. mess 'Mini yplus ' (mini yplus) ' Maxi yplus ' (maxi yplus) ;
  188.  
  189. pe=exco rv.inco.'PRESSION' 'PRES' ;
  190. pn=elno $vc (kcht $vc scal kpres pe) kpres;
  191.  
  192. Si GRAPH ;
  193. trace pn paroi ;
  194. Finsi;
  195.  
  196. mess 'Mini Pres ' (mini pn) ' Maxi Pres ' (maxi pn) ;
  197.  
  198.  
  199. FIN ;
  200.  
  201.  
  202.  
  203.  
  204.  
  205.  
  206.  

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