Télécharger la plaquette Cast3M


fr en

Télécharger fcourant.dgibi

Retour à la liste

Numérotation des lignes : 

  1. * fichier : fcourant.dgibi
  2. 'OPTI' 'ECHO' 0 ;
  3. ************************************************************************
  4. * NOM         : FCOURANT
  5. * DESCRIPTION : Teste la procédure calculant la fonction de courant
  6. *               sur un ecoulement de Poiseuille en 2D et 2D axi
  7. *
  8. *
  9. *
  10. * LANGAGE     : GIBIANE-CAST3M
  11. * AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  12. *               mél : gounand@semt2.smts.cea.fr
  13. **********************************************************************
  14. * VERSION    : v1, 21/12/2007, version initiale
  15. * HISTORIQUE : v1, 21/12/2007, création
  16. * HISTORIQUE : 2018/10/15 : test nouveaux mots-clés FCOURANT
  17. * HISTORIQUE :
  18. * HISTORIQUE :
  19. ************************************************************************
  20. * Prière de PRENDRE LE TEMPS de compléter les commentaires
  21. * en cas de modification de ce sous-programme afin de faciliter
  22. * la maintenance !
  23. ************************************************************************
  24. 'SAUTER' 2 'LIGNE' ;
  25. 'MESSAGE' ' Execution de fcourant.dgibi' ;
  26. 'SAUTER' 2 'LIGNE' ;
  27. *
  28. interact = FAUX ;
  29. graph    = FAUX ;
  30. *
  31. 'OPTION' 'DIME' 2 'ELEM' 'TRI6' ;
  32. *
  33. 'SI' ('NON' interact) ;
  34.   'OPTION' 'TRAC' 'PS' ;
  35. 'SINON' ;
  36.   'OPTION' 'TRAC' 'X' ;
  37. 'FINSI' ;  
  38. *
  39. lok = VRAI ;
  40. *
  41. * iiaxi = 1 : mode plan
  42. * iiaxi = 2 : mode axi mot-cle REGU pour FCOURANT
  43. * iiaxi = 3 : mode axi mot-cle NORM pour FCOURANT
  44. *
  45. 'REPETER' iaxi 3 ;
  46.    iiaxi = &iaxi ;
  47.    'SAUTER' 1 'LIGNE' ;
  48.    'SI' ('EGA' iiaxi 1) ;
  49.       'OPTI' 'MODE' 'PLAN' ;
  50.       'MESSAGE' 'Mode plan' ;
  51.       tdef = 'mod=PLAN' ;
  52.    'SINON' ;
  53.       'OPTI' 'MODE' 'AXIS' ;
  54.       'SI' ('EGA' iiaxi 2) ;
  55.          'MESSAGE' 'Mode axi calcul "regularise"' ;
  56.          motclefc = 'CHAI' 'REGU' ;
  57.       'SINO' ;         
  58.          'MESSAGE' 'Mode axi calcul "normal"' ;
  59.          motclefc = 'CHAI' 'NORM' ;
  60.       'FINS' ;         
  61.       tdef = 'CHAI' 'mod=AXIS mcle=' motclefc ;
  62.    'FINSI' ;
  63.    'SAUTER' 1 'LIGNE' ;
  64. *
  65. laxi = ('EGA' ('VALEUR' 'MODE') 'AXIS') ;
  66. *
  67. * Maillage
  68. *
  69. a = 1.2 ; b = 2.5 ;
  70. nx = 5 ; ny = 7 ;
  71. *nx = 10 ; ny = 14 ;
  72. dpdz = -1.5 ; mu = 0.9 ;
  73. vtoim = 0. ('*' dpdz b) ;
  74. disv = 'QUAF' ;
  75. disp = 'CENTREP1' ;
  76. *
  77. pA = 0. 0. ; pB = a 0. ; pC = a b ; pD = 0. b ;
  78. bas = 'DROIT' nx pA pB ; dro = 'DROIT' ny pB pC ;
  79. hau = 'DROIT' nx pC pD ; gau = 'DROIT' ny pD pA ;
  80. *
  81. cmt = bas 'ET' dro 'ET' hau 'ET' gau ;
  82. mt  = 'SURFACE' cmt ;
  83. _mt = 'CHANGER' mt 'QUAF' ;
  84. *
  85. $mt  = 'MODE' _mt  'NAVIER_STOKES' disv ;
  86. $hau = 'MODE' hau 'NAVIER_STOKES' disv ;
  87. *
  88. mt = 'DOMA' $mt 'MAILLAGE' ;
  89. *
  90. * Solution exacte 2D 
  91. *    u_x = 0
  92. *    u_y = vmax (1 - x^2/a^2)
  93. *          vmax = -dp/dz (a^2 / 2 mu)
  94. *    p   = dp/dz y
  95. *    psi = vmax x (1 - x^2/3a^2)
  96. *
  97. 'SI' ('NON' laxi) ;
  98.    vmax = '*' ('*' dpdz -1.) ('/' ('*' a a) ('*' mu 2)) ;
  99.    x = 'COORDONNEE' 1 mt ;
  100.    y = 'COORDONNEE' 2 ('DOMA' $mt disp) ;
  101.    cy = '+' ('*' ('/' ('*' x x) ('*' a a)) -1.) 1. ;
  102.    unx = 'NOMC' 'UY' ('*' vmax cy) ;
  103.    pnx = '*' dpdz y ;
  104.    cpsi = '+' ('*' ('/' ('*' x x) ('*' ('*' a a) 3)) -1.) 1. ;
  105.    psix = '*' vmax ('*' x cpsi) ;
  106. 'SINON' ;   
  107. *
  108. * Solution exacte 2D axi
  109. *    u_r = 0
  110. *    u_z = vmax (1 - r^2/a^2)
  111. *          vmax = -dp/dz (a^2 / 4 mu)
  112. *    p   = dp/dz y
  113. *    psi = vmax x (1 - x^2/3a^2)
  114. *
  115.    vmax = '*' ('*' dpdz -1.) ('/' ('*' a a) ('*' mu 4)) ;
  116.    r = 'COORDONNEE' 1 mt ;
  117.    z = 'COORDONNEE' 2 ('DOMA' $mt disp) ;
  118.    cz = '+' ('*' ('/' ('*' r r) ('*' a a)) -1.) 1. ;
  119. * c'est la convention EQEX, EXEC d'avoir toujours UX et UY comme noms.
  120.    unx = 'NOMC' 'UY' ('*' vmax cz) ;
  121.    pnx = '*' dpdz z ;
  122.    cpsi = '+' ('*' ('/' ('*' r r) ('*' ('*' a a) 2)) -1.) 1. ;
  123.    psix = '*' ('*' vmax PI) ('*' ('*' r r) cpsi) ;
  124. 'FINSI' ;   
  125. psix = 'NOMC' 'PSI' psix ;
  126. *
  127. * Poiseuille
  128. *
  129. rv = 'EQEX' 'ITMA' 1 'ALFA' 1. 'NITER' 1
  130.      'OPTI' 'EF' 'IMPL' 'CENTREE' disp
  131.      'ZONE' $mt  'OPER' 'NS' 1. 'UN' mu 'INCO' 'UN'
  132.      'ZONE' $mt  'OPER' 'KBBT' 1. 'INCO' 'UN' 'PN'
  133.      'ZONE' $hau 'OPER' 'TOIM' vtoim 'INCO' 'UN'
  134.      'CLIM' 'UN' 'UIMP' dro 0. 'UN' 'VIMP' dro 0.
  135.      'CLIM' 'UN' 'UIMP' gau 0. 
  136.    ;
  137. rv . 'INCO' = 'TABLE' 'INCO' ;
  138. RV . 'INCO' . 'UN'   = 'KCHT' $mt 'VECT' 'SOMMET' (0. 0.) ;
  139. RV . 'INCO' . 'PN' = 'KCHT' $mt 'SCAL' disp 0. ;
  140. *
  141. EXEC rv ;
  142. *
  143. un  = rv . 'INCO' . 'UN' ;
  144. pn  = rv . 'INCO' . 'PN' ;
  145. *clpsi = 'MANUEL' 'CHPO' pA 1 'SCAL' 0. ;
  146. blpsi = 'BLOQ' 'T' pA ;
  147. 'SI' laxi ;
  148.    psi = FCOURANT mt un blpsi motclefc ;
  149. 'SINO' ;
  150.    psi = FCOURANT mt un blpsi ;
  151. 'FINS' ;
  152. * 
  153. 'SI' graph ;
  154.    titg = 'CHAINE' ' ' tdef ;
  155.    vun  = 'VECT' un 'UX' 'UY' 'JAUN' ;
  156.    tit = 'CHAINE' 'Vitesse' titg ;
  157.    'TRACER' vun mt cmt 'TITR' tit ;
  158.    ppn = 'ELNO' $mt pn disp ;
  159.    tit = 'CHAINE' 'Pression' titg ;
  160.    'TRACER' ppn mt  'TITR' tit ;
  161.    tit = 'CHAINE' 'Fonction de courant' titg ;
  162.    'TRACER' psi mt  'TITR' tit ;
  163. 'FINSI' ;
  164. *
  165. erru = 'MAXIMUM' ('-' un unx) 'ABS' ;
  166. errp = 'MAXIMUM' ('-' pn pnx) 'ABS' ;
  167. errf = 'MAXIMUM' ('-' psi psix) 'ABS' ;
  168. *
  169. 'MESSAGE' ('CHAINE' 'erru = ' erru) ;
  170. 'MESSAGE' ('CHAINE' 'errp = ' errp) ;
  171. 'MESSAGE' ('CHAINE' 'errf = ' errf) ;
  172. *
  173. * le 21/12/2007 sur ma machine  16/10/2018 
  174. * Mode plan
  175. *   erru =  4.21885E-15   erru =  3.33067E-15
  176. *   errp =  3.15303E-14   errp =  2.08722E-14
  177. *   errf =  6.36767E-04   errf =  6.36767E-04
  178. * Mode axi regularise
  179. *   erru =  2.66454E-15   erru =  4.28546E-14
  180. *   errp =  2.62013E-14   errp =  1.14575E-13
  181. *   errf =  5.45703E-03   errf =  5.45703E-03
  182. * Mode axi normal
  183. *                         erru =  4.28546E-14
  184. *                         errp =  1.14575E-13
  185. *                         errf =  5.89719E-03
  186. *
  187. tolup = 1.D-12 ;
  188. 'SI' ('EGA' iiaxi 1) ;
  189.    tolf = 7.D-4 ;
  190. 'FINS' ;
  191. 'SI' ('EGA' iiaxi 2) ;
  192.    tolf = 6.D-3 ;
  193. 'FINSI' ;
  194. 'SI' ('EGA' iiaxi 3) ;
  195.    tolf = 6.D-3 ;
  196. 'FINSI' ;
  197. *
  198. tst = ('<' erru tolup) ;
  199. 'SI' ('NON' tst) ;
  200.    chmes = 'CHAINE' '!!! Erreur, on na pas erru < ' tolup ;
  201.    'MESSAGE' chmes ;
  202. 'FINSI' ;
  203. lok = lok 'ET' tst ;   
  204. *
  205. tst = ('<' errp tolup) ;
  206. 'SI' ('NON' tst) ;
  207.    chmes = 'CHAINE' '!!! Erreur, on na pas errp < ' tolup ;
  208.    'MESSAGE' chmes ;
  209. 'FINSI' ;
  210. lok = lok 'ET' tst ;   
  211. *
  212. tst = ('<' errf tolf) ;
  213. 'SI' ('NON' tst) ;
  214.    chmes = 'CHAINE' '!!! Erreur, on na pas errf < ' tolf ;
  215.    'MESSAGE' chmes ;
  216. 'FINSI' ;
  217. lok = lok 'ET' tst ;
  218. *
  219. 'FIN' iaxi ;
  220. *
  221. 'SAUTER' 2 'LIGNE' ;
  222. 'SI' lok ;
  223.    'MESSAGE' 'Tout sest bien passe' ;
  224. 'SINON' ;
  225.    'MESSAGE' 'Il y a eu des erreurs' ;
  226. 'FINSI' ;   
  227. 'SAUTER' 2 'LIGNE' ;
  228. *
  229. 'SI' interact ;
  230.    'OPTION' 'DONN' 5 'ECHO' 1 ;
  231. 'FINSI' ;   
  232. 'SI' ('NON' lok) ;
  233.    'ERREUR' 5 ;
  234. 'FINSI' ;   
  235. *
  236. * End of dgibi file FCOURANT
  237. *
  238. 'FIN' ;
  239.  
  240.  
  241.  
  242.  
  243.  
  244.  

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