Télécharger oscicyl2.dgibi

Retour à la liste

Numérotation des lignes :

  1. ************************************************************************
  2. * NOM : OSCICYL2
  3. * DESCRIPTION : Cylindre oscillant dans une cavité fermée
  4. *
  5. * Deux paramètres : le nombre de Reynolds (convection/viscosité)
  6. * le nombre de Keulegan-Carpenter (convection/masse,trainee)
  7. *
  8. *
  9. * Bibliographie : thèse de Marion Duclercq
  10. * avec Daniel BROC
  11. *
  12. * LANGAGE : GIBIANE-CAST3M
  13. * AUTEUR : Stephane GOUNAND et Dihya MEZI (CEA/DES/ISAS/DM2S/SEMT/LTA)
  14. * mel : stephane.gounand@cea.fr
  15. **********************************************************************
  16. * VERSION : v1, 28/01/2022, version initiale
  17. * HISTORIQUE : v1, 28/01/2022, création
  18. * HISTORIQUE :
  19. * HISTORIQUE :
  20. ************************************************************************
  21. *
  22. *
  23. 'OPTION' 'DIME' 2 'ELEM' 'QUA8' 'ECHO' 0 ;
  24. *'OPTI' 'ISOV' suli ;
  25.  
  26. graph = faux ;
  27. * condition vitesse normale
  28. l_cl = vrai ;
  29. ****Definitions des parametres des longueurs
  30.  
  31. Rayon = 20.E-2;
  32. Long = 13 * Rayon ;
  33. Haut = 10 * Rayon ;
  34.  
  35. * Parametres pour le maillage
  36.  
  37. nx=26;
  38. ny=14;
  39. nc=6;
  40. *nx = nx * 2 ;
  41. *ny = ny * 2 ;
  42. *nc = nc * 2 ;
  43.  
  44. ********** Creation du contour exterieur
  45. pA= 0. 0. ; pB=Long 0. ; pC= Long Haut ; pD= 0. Haut ;
  46. bas = 'DROI' nx pA pB;
  47. dro= 'DROI' ny pB pC;
  48. hau= 'DROI' nx pC pD;
  49. gau = 'DROI' ny pD pA;
  50. **** Creation du cercle
  51. pcen= (Long * 0.5) (Haut * 0.5) ;
  52. pE = ((0.5 * Long) '-' Rayon) (0.5 * Haut);
  53. pF = ((0.5 * Long) ) ((0.5 * Haut) '+'Rayon );
  54. pG= ((0.5 * Long)'+' Rayon ) ((0.5 * Haut) );
  55. pH=((0.5 * Long) ) ((0.5 * Haut) '-' Rayon );
  56. CE1 = 'CERCLE' nc pE pcen pF;
  57. CE2= 'CERCLE' nc CE1 pcen pG;
  58. CE3= 'CERCLE' nc CE2 pcen pH;
  59. CE= 'CERCLE' nc CE3 pcen pE;
  60. * Contour du maillage exterieur
  61. cmt1= bas 'ET' dro 'ET' hau 'ET' gau ;
  62. *** Contour total avec cercle
  63. cmt= bas 'ET' dro 'ET' hau 'ET' gau 'DIFF' CE;
  64. mt = 'SURFACE' cmt;
  65. 'ELIMINATION' mt 1.E-9;
  66. * Changer le maillage en un maillage quadratique fluide
  67. _mt='CHAN' mt 'QUAF';
  68. *'TRAC' _mt nclk ;
  69. ************************************************************************
  70. * champ des normales CE
  71. CE2 = changer CE seg2 ;
  72. mod_cyl2 = mode CE2 mecanique elastique coq2 ;
  73. norm_cyl = pression mod_cyl2 1. coque normale ;
  74. norm_cyl = norm_cyl /
  75. ((psca norm_cyl norm_cyl (mots fx fy) (mots fx fy)) ** 0.5) ;
  76. v_norm = vecteur norm_cyl 1.e-2 fx fy jaune ;
  77. si graph ;
  78. titre (et tit_ref ' / champ des normales ') ;
  79. trac v_norm CE nclk;
  80. finsi ;
  81. *
  82. n_CE = NOMC ('ENLE' norm_cyl 'MZ') (mots 'FX' 'FY') (mots 'UX' 'UY') ;
  83. mn_CE = extr n_CE mail ;
  84. ie1 = 0 ;
  85. repeter b_blCE (nbel mn_CE) ;
  86. ie1 = ie1 + 1 ;
  87. poie1 = 'POINT' mn_CE ie1 ;
  88. ch_u = 'REDU' n_CE poie1 ;
  89. mbl_u = 'RELA' ch_u ;
  90. fbl_u = 'DEPI' mbl_u ('EXTR' ch_u 'UX' poie1) ;
  91. si (ie1 ega 1) ;
  92. mbl_CE = mbl_u ;
  93. fbl_CE = fbl_u ;
  94. sinon ;
  95. mbl_CE = mbl_u et mbl_CE ;
  96. fbl_CE = fbl_u et fbl_CE ;
  97. finsi ;
  98. fin b_blCE ;
  99. fbl_CEv = fbl_CE ;
  100.  
  101. **** Conditions de blocage
  102. *** Vitesse
  103. *mblox = 'BLOQUE' 'UX' cmt1 ;
  104. mblox = 'BLOQUE' 'UX' (dro et gau) ;
  105. *mbloy = 'BLOQUE' 'UY' cmt1 ;
  106. mbloy = 'BLOQUE' 'UY' (bas et hau) ;
  107. mblocx = 'BLOQUE' 'UX' CE ;
  108. mblocy = 'BLOQUE' 'UY' CE ;
  109.  
  110. **** Physique
  111. *** nombre de Reynolds de 100
  112. *** nombre de Keulegan-Carpenter de 2
  113. Re = 100. ; rho1 = 1.e3 ; mu1 = 1.e0 ;
  114. nre= 1 '/' Re ;
  115. KC = 2. ;
  116.  
  117. **** Numérique
  118. *** Pas de temps (cf. thèse Marion Duclercq)
  119. 'SI' ('>' KC 2.) ;
  120. dt = 0.05 '/' KC ;
  121. 'SINO' ;
  122. dt = 0.025 ;
  123. 'FINS' ;
  124. *dt = dt * 4. ;
  125.  
  126.  
  127. mns1 = mode _mt 'NAVIER_STOKES' 'NLIN' 'QUAF' 'LINE' 'CHPO' 'GAP7' ;
  128. cns1 = mate mns1 RHO rho1 MU mu1 ;
  129. cns2 = cns1 et (mate mns1 FREQ (dt**(-1))) ;
  130. *** Pression
  131. mapcen = extr mns1 'CENT' ;
  132. pp2 = 'POIN' mapcen 'PROC' (0. (Haut/2)) ;
  133. mblop = 'BLOQ' 'LX1' pp2 ;
  134.  
  135.  
  136.  
  137. mtbid = changer _mt line ;
  138. mbid = mode mtbid mecanique elastique ;
  139. cbid = mate mbid youn 1.e10 nu 0.3 rho 2.;
  140.  
  141. **** Matrices fixes du probleme
  142. *** Rigidité (objet RIGIDITE)
  143. * mrigi = 'NLIN' gdisc _mt A B mot1 mot2 ;
  144. nrigi = rigi (mns1 et mbid) (cns1 et cbid) ;
  145. rigbid = rigi (mbid) (cbid) ;
  146. *** Contrainte div U
  147. ndiv = rigi 'DIVE' (mns1 et mbid) (cns1 et cbid) ;
  148. *** Masse pondérée par le pas de temps
  149. nmass = mass (mns1 et mbid) (cns2 et cbid) ;
  150. ************************************************************************
  151.  
  152.  
  153. **** Conditions initiales
  154. U0 = 'MANUEL' 'CHPO' _mt 2 'UX' 0. 'UY' 0. ;
  155. UNM = U0 ;
  156. UN = U0 ;
  157. solu= UN;
  158. EPS= 1.0D-9;
  159. tps= 0. ;
  160.  
  161. *
  162. * mouvement impose au cylindre
  163. abs_t = prog 0. pas Dt (20*Dt) ;
  164. *abs_t = prog 0. pas Dt (70*Dt) ;
  165. ord_t = 'SIN' (360. '*' abs_t) * Rayon * 1.e-1 ;
  166. ev_t = evol vert manu ' t' abs_t ord_t ;
  167. absx_CE = ((coor 1 CE) - (coor 1 pcen))/Rayon ; absc_CE = abs absx_CE ;
  168. absy_CE = ((coor 2 CE) - (coor 2 pcen))/Rayon ;
  169. si l_cl ;
  170. fbl_CE = fbl_CEv ;
  171. sinon ;
  172. fbl_CE = 'DEPI' mblocx ('NOMC' 'UX' absc_CE) ;
  173. finsi ;
  174. char_t = 'CHAR' 'DIMP' fbl_CE ev_t ;
  175. unmax = maxi fbl_CE ;
  176.  
  177. *
  178. t_fltub = table 'PASAPAS' ;
  179. t_fltub . modele = mns1 ;
  180. t_fltub . caracteristiques = cns1 ;
  181. si l_cl ;
  182. t_fltub . blocages_mecaniques = mblox 'ET'
  183. mbloy 'ET' mbl_CE et mblop ;
  184.  
  185. sinon ;
  186. t_fltub . blocages_mecaniques = mblox 'ET'
  187. mbloy 'ET' mblop 'ET' mblocx 'ET' mblocy ;
  188. finsi ;
  189. t_fltub . chargement = char_t ;
  190. t_fltub . temps_calcules = abs_t ;
  191. *t_fltub . temps_calcules = prog 0. Dt ;
  192.  
  193. pasapas t_fltub ;
  194.  
  195. ************************************************************************
  196. graph = faux ;
  197. nomvit = 'EXTR' mns1 'DEPL' ;
  198. t_pap = t_fltub ;
  199. dep_0 = t_pap . vitesses_fluide . 1 * 0. ;
  200. abs_t = prog ; of_x = prog ; of_y = prog ;
  201. op_x = prog ; op_y = prog ;
  202. f0_x = 0. ; f0_y = 0. ; p0_max = 0. ; p0_min = 0. ;
  203. ux_CE = 0. ;
  204.  
  205. repeter BTRA ((dime t_pap . temps) - 1) ;
  206. solu = t_pap . vitesses_fluide . &btra ;
  207. unpip = 'EXCO' solu nomvit ;
  208. pnp = t_pap . pression_fluide . &btra ;
  209. vun = 'VECTEUR' (unpip/unmax/3) nomvit 'ORAN' ;
  210.  
  211. * trainee et portance defini par les reactions
  212. si l_cl ;
  213. ch_pr0 = 'REAC' mbl_CE solu ;
  214. * maillage CE régulier
  215. f0_CE = 'RESU' ch_pr0 ;
  216. f_x = 'EXTR' f0_CE 'FX' ('POINT' 1 ('EXTR' f0_CE 'MAIL')) ;
  217. f_y = 'EXTR' f0_CE 'FY' ('POINT' 1 ('EXTR' f0_CE 'MAIL')) ;
  218. sinon ;
  219. f_x = 'MAXI' ('RESU' ('REAC' mblocx solu)) ;
  220. f_y = 'MAXI' ('RESU' ('REAC' mblocy solu)) ;
  221. finsi ;
  222. f_x = f_x * (-1) ; f_y = f_y * (-1) ;
  223. abs_t = abs_t et t_pap . temps . &btra ;
  224. of_x = of_x et (prog f_x) ; of_y = of_y et (prog f_y) ;
  225. f0_x = maxi (prog (abs f_x) f0_x) ; f0_y = maxi (prog (abs f_y) f0_y) ;
  226.  
  227. * resultante pression
  228. pre_CE = redu pnp CE ;
  229. p_max = maxi pre_CE ; p0_max = maxi (prog p_max p0_max) ;
  230. p_min = mini pre_CE ; p0_min = mini (prog p_min p0_min) ;
  231. pre_x = pre_cE * absx_CE ; pre_y = pre_CE * absy_CE ;
  232. fp_x = intg (evol chpo CE pre_x) ; fp_y = intg (evol chpo CE pre_y) ;
  233. op_x = op_x et (prog fp_x) ; op_y = op_y et (prog fp_y) ;
  234.  
  235. fin BTRA ;
  236.  
  237. of_x = of_x / f0_x ; of_y = of_y / f0_y ;
  238.  
  239. si (graph) ;
  240.  
  241. repeter BTRA ((dime t_pap . temps) - 1) ;
  242. solu = t_pap . vitesses_fluide . &btra ;
  243. unpip = 'EXCO' solu nomvit ;
  244. pnp = t_pap . pression_fluide . &btra ;
  245. vun = 'VECTEUR' (unpip/unmax/3) nomvit 'ORAN' ;
  246. chf_cer = manu chpo pcen 2 'FX' (extr of_x &btra) 'FY' (extr of_y &btra) ;
  247. vf_cer = vecteur 0.7 chf_cer 'FX' 'FY' 'TURQ' ;
  248.  
  249. v_visu = vun et vf_cer ; ma_visu = _mt et pcen ;
  250.  
  251. titres 'pression / vitesses temps = ' t_pap . temps . &btra ;
  252. trac pnp v_visu ma_visu (contour mt) nclk;
  253.  
  254. 'SI' ('EGA' ('TYPE' mz_CE) 'MAILLAGE') ;
  255. unz_CE = redu unpip mz_CE ;
  256. vunz_CE = 'VECTEUR' (unz_CE/unmax) nomvit 'ORAN' ;
  257. trac pnp vunz_CE mz_CE(contour mz_CE) ;
  258. 'FINSI' ;
  259.  
  260. si (ega &btra 1) ; def_1 = defo ma_visu dep_0 v_visu pnp ;
  261. sinon ;
  262. def_1 = def_1 et (defo ma_visu dep_0 v_visu pnp) ;
  263. finsi ;
  264.  
  265.  
  266. fin BTRA ;
  267.  
  268. evf_x = 'COUL' vert ('EVOL' 'MANU' abs_t of_x) ;
  269. evf_y = 'COUL' bleu ('EVOL' 'MANU' abs_t of_y) ;
  270. 'DESS' (evf_x et evf_y) 'TITR' 'f_x(t) parallele vert / f_y(t) orthogonal bleu (normes)' ;
  271. mess 'reactions x maxi : ' f0_x ' y maxi : ' f0_y ;
  272.  
  273. evp_x = 'COUL' rouge ('EVOL' 'MANU' abs_t op_x) ;
  274. evp_y = 'COUL' jaune ('EVOL' 'MANU' abs_t op_y) ;
  275. 'DESS' ((evf_x * f0_x) et evp_x) 'TITR' 'reaction vert / pression rouge en x ' ;
  276. 'DESS' ((evf_y * f0_y) et evp_y) 'TITR' 'reaction bleu / pression jaune en y ' ;
  277.  
  278.  
  279. d_cop = (p0_max - p0_min) / 16 ; l_cop = prog p0_min pas d_cop p0_max ;
  280. titre 'vitesses / pression / effort cylindre dt = ' dt ' t max = ' (maxi abs_t) ;
  281. trac anim def_1 l_cop ;
  282. *
  283.  
  284. finsi ;
  285.  
  286. ef_x = maxi (abs ((of_x*f0_x) - op_x)) ; mess 'ecart ' ef_x ;
  287. si (ef_x < 1.e-3) ; erre 0 ; sinon ; erre 5 ; finsi ;
  288.  
  289.  
  290.  
  291. *opti donn 5 ;
  292. *
  293. * End of dgibi file oscillations
  294. *
  295. 'FIN' ;
  296.  
  297.  
  298.  
  299.  

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