Télécharger vibr12.dgibi

Retour à la liste

Numérotation des lignes :

  1. * fichier : vibr12.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4.  
  5. ************************************************************************
  6. * *
  7. * Mots-cles : Vibrations, calcul modal, modes complexes, *
  8. * interaction fluide-structure, instabilite *
  9. * *
  10. * TEST : VIBR12 *
  11. * *
  12. * CALCUL D INSTABILITE FLUIDE-ELASTIQUE SOUS ECOULEMENT AXIAL CONFINE *
  13. * 2 CONFIGURATIONS : ENCASTRE-LIBRE ou AA *
  14. * *
  15. * Ref : GIBERT Page 409-410 *
  16. * GREGORY, PAIDOUSSIS, 1966, "Unstable oscillation of tubular *
  17. * cantilevers conveying fluid " *
  18. * *
  19. * ^ Z *
  20. * -> | *
  21. * -> V +--------------------+--------------------+ --> X *
  22. * -> P0 P2 *
  23. * *
  24. * *
  25. ************************************************************************
  26.  
  27. ************************************
  28. * OPTIONS *
  29. ************************************
  30. OPTI DIME 3 ELEM SEG2 ;
  31. OPTI EPSI LINEAIRE ;
  32.  
  33. * Conditions limites : 'EL' ou 'AA'
  34. CL = mot 'AA' ;
  35. CL = mot 'EL' ;
  36. prefix = CHAI 'vibr12-' CL;
  37.  
  38. * traces
  39. * GRAPH = VRAI ;
  40. * OPTI TRAC PSC EPTR 10 POTR 'HELVETICA_16' 'FTRA' (CHAI prefix'.ps');
  41. GRAPH = FAUX ;
  42. OPTI TRAC PSC EPTR 10 POTR 'HELVETICA_16' 'FTRA' (CHAI prefix'.ps');
  43.  
  44.  
  45.  
  46. ************************************
  47. * Caracteristiques du systeme *
  48. ************************************
  49.  
  50. *--- Materiau Acier
  51. Rho_A = 8100. ;
  52. You_A = 2.0e11 ;
  53. Nu_A = 0.3 ;
  54.  
  55. *--- Longueurs : Longueur Tube = 2L
  56. L = .5 ;
  57.  
  58. *--- Diametre int. Tube
  59. Di = 10.0e-3 ;
  60.  
  61. *--- Rapport Beta = Mf / (Mf+ms)
  62. Beta = 0.2 ;
  63.  
  64. *--- choix d'une methode de calcul pour :
  65. * | Rhof : Masse volumique fluide
  66. * | De : diametre externe
  67. * sachant que : | ms = Mf * (1-Beta) / Beta
  68. * | (De^2-Di^2) = Rhof*Di^2/Rhos * (1-Beta)/Beta
  69. Meth = 1 ;
  70. SI (EGA Meth 1) ;
  71. * On calcul le Diametre externe De pour respecter le rapport Beta
  72. * a masse volumique fluide fixe
  73. Rhof = 1000. ;
  74. De2 = ( (((Rhof/Rho_A) * ( (1-Beta) / Beta ) )) + 1 ) * (Di*Di) ;
  75. De = De2**0.5 ;
  76. SINON ;
  77. * On calcul Rhof pour respecter Beta ;
  78. * a diametre exterieur fixe
  79. De = 12e-3 ;
  80. Rhof = (((( (De*De) - (Di*Di)) / ((1-Beta)/Beta) ) ) / (Di*Di) ) *Rho_A ;
  81. FINSI ;
  82.  
  83. *--- Epaisseur du tube
  84. Ep = (De - Di) / 2 ;
  85.  
  86. Mess 'Di=' Di ' De=' De 'Rhof=' Rhof ;
  87.  
  88. *--- Section Tube
  89. S_T = (pi/4)*((De**2)-(Di**2)) ;
  90.  
  91. *--- Moments quadratiques
  92. Ie_T= (pi/64)*((De**4)-(Di**4)) ;
  93.  
  94. *--- Masse ajoutee lineique interne
  95. S_I = (pi/4.0)*(Di**2) ;
  96. Me = Rhof*S_I ;
  97. Rho_eq = Me/S_T ;
  98.  
  99. * Nbre de modes reels et amortissement modaux en % (cf operateur AMOR)
  100. Nbmod = 4 ;
  101. L_AMOR = PROG Nbmod * 0.0 ;
  102.  
  103. * Facteurs pour adimensionner [Gibert,p.408] :
  104. * V : VR = RVR*V
  105. * w : OMEG = Rw*w
  106. RVR = (2.0*L) * (( ((Rho_A*S_T) + Me) / (You_A*Ie_T) )**0.5) ;
  107. Rw = (2.0*L)* RVR ;
  108. * Autre adimensionnement [Axisa,Paidoussis]
  109. si (ega CL 'EL');
  110. RVR = (2.0*L) * (( Me / (You_A*Ie_T) )**0.5) ;
  111. finsi;
  112.  
  113. * Juste pour verif
  114. Beta_c= Rhof*S_I / ((Rhof*S_I) + (Rho_A*S_T)) ;
  115.  
  116. *--- Vitesse axial du fluide
  117. pas_VR = .1 ;
  118. VRmax = 12. ;
  119. VRmin = 0. ;
  120. LVR = PROG VRmin PAS pas_VR VRmax ;
  121. LVa = (1.0 * LVR ) / RVR ;
  122.  
  123. mess ' Rho A Rho_eq Beta Beta_c Rw RVr De' ;
  124. mess ' ' Rho_A ' ' Rho_eq ' ' Beta ' ' Beta_c ' ' Rw ' ' RVR ' ' De;
  125.  
  126.  
  127. **************************************
  128. * Maillage du systeme *
  129. **************************************
  130.  
  131. P0 = (-1.0*L) 0. 0. ;
  132. P2 = P0 PLUS ((2.0*L) 0. 0.) ;
  133.  
  134. * Lfin : longueur objectif d'1 element
  135. * L1 : longueur reelle
  136. *Lfin = 0.001 ;
  137. Lfin = 0.0025 ;
  138. nL = ENTI 'PROCHE' (2.0*L / Lfin);
  139. Tube = P0 DROI nL P2 ;
  140. * Longueur reelle des elements
  141. * L1 = (MESURE Tube) / nL ;
  142. L1 = 2.0*L / nL;
  143.  
  144. **************************************
  145. * Modele et Materiau *
  146. **************************************
  147.  
  148. Modl1 = MODE Tube MECANIQUE ELASTIQUE ISOTROPE 'TUYA' ;
  149. Matr1 = MATE Modl1 'YOUN' You_A 'NU' Nu_A 'RHO' (Rho_A + Rho_eq) 'EPAI' Ep 'RAYO' (De/2);
  150.  
  151.  
  152. **************************************
  153. * Masses, Rigidites ... du systeme *
  154. **************************************
  155.  
  156. Mstru = MASS Modl1 Matr1 ;
  157. Kstru = RIGI Modl1 Matr1 ;
  158.  
  159. *--- Conditions aux limites : Probleme plan EL ou EA
  160. Cl1 = BLOQ 'UX' Tube ;
  161. Cl2 = BLOQ 'UY' Tube ;
  162. Cl3 = BLOQ 'RX' Tube ;
  163. Cl4 = BLOQ 'RZ' Tube ;
  164. * EL : Encastre-Libre
  165. SI (EGA CL 'EL') ;
  166. Cl5 = BLOQ 'RY' P0 ;
  167. SINON ;
  168. * AA : Appuye-Appuye
  169. SI (EGA CL 'AA') ;
  170. Cl5 = BLOQ 'UZ' P2 ;
  171. SINON ;
  172. MESS 'CL doit etre AA ou EL !';
  173. ERRE 21;
  174. FINSI ;
  175. FINSI ;
  176. Cl6 = BLOQ 'UZ' P0 ;
  177. Kstru = Kstru ET Cl1 ET Cl2 ET Cl3 ET Cl4 ET CL5 ET CL6 ;
  178.  
  179.  
  180. *****************************************
  181. * Analyse modale *
  182. *****************************************
  183.  
  184. *Basr = VIBR 'PROC' (PROG 1.) (LECT Nbmod) Kstru Mstru ;
  185. Basr = VIBR 'IRAM' 1. Nbmod Kstru Mstru ;
  186. Basr.'MODES' = NNOR (Basr.'MODES') (MOTS 'UZ') ;
  187. Basrm = Basr.'MODES' ;
  188.  
  189. OPTI OEIL (0.0 -1000.0 0.0) ;
  190. POSTVIBR Basr ;
  191.  
  192.  
  193. ***********************************
  194. * Rigidites suplementaires *
  195. ***********************************
  196.  
  197. PAfu1 = (PROG -.5 0. 0. (L1/10.) .5 0. 0. (-1*L1/10.)) ;
  198. PAfu2 = (PROG 0. -.5 (-1*L1/10.) 0. 0. .5 (L1/10.) 0.) ;
  199. PAfu3 = (PROG 0. (L1/10.) 0. 0. 0. (-1*L1/10.) (-1*L1*L1/60.) 0.) ;
  200. PAfu4 = (PROG (-1*L1/10.) 0. 0. 0. (L1/10.) 0. 0. (-1*L1*L1/60.)) ;
  201. PAfu5 = (PROG -.5 0. 0. (-1*L1/10.) .5 0. 0. (L1/10.)) ;
  202. PAfu6 = (PROG 0. -.5 (L1/10.) 0. 0. .5 (-1*L1/10.) 0.) ;
  203. PAfu7 = (PROG 0. (-1*L1/10.) (L1*L1/60.) 0. 0. (L1/10.) 0. 0.) ;
  204. PAfu8 = (PROG (L1/10.) 0. 0. (L1*L1/60.) (-1*L1/10.) 0. 0. 0.) ;
  205. PAfu = PAfu1 ET PAfu2 ET PAfu3 ET PAfu4 ET PAfu5 ET PAfu6 ET PAfu7 ET PAfu8 ;
  206. Afu = MANU 'RIGI' 'TYPE' 'AMORTISSEMENT'
  207. Tube (MOTS 'UY' 'UZ' 'RY' 'RZ') 'QUEL' PAfu ;
  208.  
  209. PKfu1 = (PROG (-6./(5.*L1)) 0. 0. -1.1 (6./(5.*L1)) 0. 0. -.1) ;
  210. PKfu2 = (PROG 0. (-6./(5.*L1)) 1.1 0. 0. (6./(5.*L1)) .1 0.) ;
  211. PKfu3 = (PROG 0. .1 (-2.*L1/15.) 0. 0. -.1 (L1/30.) 0.) ;
  212. PKfu4 = (PROG -.1 0. 0. (-2.*L1/15.) .1 0. 0. (L1/30.)) ;
  213. PKfu5 = (PROG (6./(5.*L1)) 0. 0. .1 (-6./(5.*L1)) 0. 0. 1.1) ;
  214. PKfu6 = (PROG 0. (6./(5.*L1)) -.1 0. 0. (-6./(5.*L1)) -1.1 0.) ;
  215. PKfu7 = (PROG 0. .1 (L1/30.) 0. 0. -.1 (-2.*L1/15.) 0.) ;
  216. PKfu8 = (PROG -.1 0. 0. (L1/30.) .1 0. 0. (-2.*L1/15.)) ;
  217. PKfu = PKfu1 ET PKfu2 ET PKfu3 ET PKfu4 ET PKfu5 ET PKfu6 ET PKfu7 ET PKfu8 ;
  218. Kfu = MANU 'RIGI' 'TYPE' 'RIGIDITE'
  219. Tube (MOTS 'UY' 'UZ' 'RY' 'RZ') 'QUEL' PKfu ;
  220.  
  221.  
  222. list 'RESUME' Afu;
  223. list 'RESUME' Kfu;
  224.  
  225.  
  226. *****************************************
  227. * Mise en forme du tableau de resultats *
  228. *****************************************
  229.  
  230. Mod = TABLE ;
  231. Mod.Reel = TABLE ;
  232. Mod.Imag = TABLE ;
  233. Repeter Bloc6 (2*Nbmod) ;
  234. i = &Bloc6 ;
  235. Mod.Reel.i = PROG ;
  236. Mod.Imag.i = PROG ;
  237. Fin Bloc6 ;
  238.  
  239. *****************************************
  240. * Projection des matrices *
  241. *****************************************
  242.  
  243. *--- Struture Mt,Ksp
  244. Mt = PJBA MStru Basr ;
  245. Ksp = PJBA Kstru Basr ;
  246. *--- Amortissement modale Ct
  247. Ct = AMOR Basr L_AMOR ;
  248.  
  249. *--- Matrice amortissement due a l ecoulement
  250. * a multiplier par 2MV pour Afp et MV² pour Kfp
  251. Afp = PJBA Afu Basr ;
  252. Kfp = PJBA Kfu Basr ;
  253.  
  254. NbModC = 2*Nbmod;
  255. mocomp = MOTS 'ALFA';
  256. numodC = LECT 1 'PAS' 1 NbModC;
  257.  
  258. * ===================================== *
  259. * Boucle sur la vitesse axiale *
  260. * ===================================== *
  261. REPE b_Vax (DIME LVa) ;
  262. k = &b_Vax ;
  263. Va= Extr LVa k ;
  264. Vr= Extr LVR k ;
  265.  
  266. Mess ' ' ;
  267. Mess '#' k ': Vitesse axiale =' Va ' Vitesse reduite =' Vr ;
  268. Mess '-------------------------------------------------------------------------------' ;
  269. Mess ' Mode fn_d Amort. P. Reeln1 P. Imag. Stabilite ' ;
  270.  
  271. * Matrices fluides totale : Amortissement et Raideur
  272. At = Ct ET ((2.*Me*Va)*Afp) ;
  273. Kt = Ksp ET (Me*Va*Va*Kfp);
  274.  
  275. * Calcul des modes complexes
  276. Basc = VIBC Mt Kt At ;
  277. Bascm = Basc . 'MODES' ;
  278.  
  279. *--- Extraction des frequences
  280. Prf = PROG ;
  281. Pif = PROG ;
  282. REPE Bloc2 NbModC ;
  283. i = &Bloc2 ;
  284. Modn = Bascm.i ;
  285. Prf = Prf ET (PROG Modn.'FREQUENCE_REELLE') ;
  286. Pif = Pif ET (PROG Modn.'FREQUENCE_IMAGINAIRE') ;
  287. FIN bloc2 ;
  288.  
  289. *--- Tri croissant des frequences
  290. * NbModC = DIME Prf ;
  291. * Prfc Pifc numodC = ORDO Prf Pif numodC0;
  292. * ou plutot suivi des modes (auteur : Benoit Prabel)
  293. si (ega k 1);
  294. * 1er pas : on ordonne par fr croissant
  295. Prfc Pifc numodCc = ORDO Prf Pif numodC;
  296. Bascmc = COPI Bascm;
  297. repe b0 NbModC;
  298. Bascmc . &b0 = Bascm . (extr numodCc &b0);
  299. fin b0;
  300. sinon;
  301. xcou01 = prog;
  302. * boucle sur les modes du pas precedent --------------------
  303. repe b0 NbModC;
  304. DEFR0 = Bascm0 . &b0 . 'DEFORMEE_MODALE_REELLE';
  305. DEFI0 = Bascm0 . &b0 . 'DEFORMEE_MODALE_IMAGINAIRE';
  306. frqR0 = Bascm0 . &b0 . 'FREQUENCE_REELLE';
  307. frqI0 = Bascm0 . &b0 . 'FREQUENCE_IMAGINAIRE';
  308. * boucle sur les modes du pas actuel ----------------
  309. repe b1 NbModC;
  310. DEFR1 = Bascm . &b1 . 'DEFORMEE_MODALE_REELLE';
  311. DEFI1 = Bascm . &b1 . 'DEFORMEE_MODALE_IMAGINAIRE';
  312. frqR1 = Bascm . &b1 . 'FREQUENCE_REELLE';
  313. frqI1 = Bascm . &b1 . 'FREQUENCE_IMAGINAIRE';
  314. * ecart en frequence^2
  315. dfrq1 = ((frqR1 - frqR0)**2) + ((frqI1 - frqI0)**2);
  316. ref1 = ((frqR1**2) + (frqI1**2)) + ((frqR0**2) + (frqI0**2));
  317. * produit scalaire des vecteurs propres
  318. cos1 = (XTY DEFR1 DEFR0 mocomp mocomp)
  319. + (XTY DEFI1 DEFI0 mocomp mocomp);
  320. sin1 = (XTY DEFI1 DEFR0 mocomp mocomp)
  321. - (XTY DEFR1 DEFI0 mocomp mocomp);
  322. modul1 = (cos1**2) + (sin1**2);
  323. xsuiv1 = (abs (1. - modul1)) + (0.5*dfrq1 / ref1);
  324. xcou01 = xcou01 et xsuiv1;
  325. fin b1;
  326. * fin de boucle sur les modes du pas actuel ---------
  327. fin b0;
  328. * fin de boucle sur les modes du pas precedent -------------
  329. * SUIVI = apprairage
  330. * = minimisation du cout (cad des distances) par permutation
  331. couttot numodCc Prfc Pifc
  332. = ORDO numodC Prf Pif 'COUT' xcou01 ;
  333. * on modifie la base modale
  334. Bascmc = COPI Bascm;
  335. repe b0 NbModC;
  336. Bascmc . &b0 = Bascm . (extr numodCc &b0);
  337. fin b0;
  338. finsi;
  339. * pour le prochain pas
  340. Bascm0 = Bascmc;
  341.  
  342. *--- Sauvegarde
  343. REPE Bloc5 NbModC ; i = &Bloc5 ;
  344. f = EXTR Prfc i ;
  345. if = EXTR Pifc i ;
  346.  
  347. Msg = 'STAB.' ;
  348. SI ((< if -1.D-10) ET (>EG f 0.) ET (f < 1e-5 )) ;
  349. Msg = 'INST.' ;
  350. SI (EGA VSauv 'N') ;
  351. Vcrit = Va_CO - (p_fin/2.) ;
  352. VSauv = 'O' ;
  353. MESS 'bipbip' ;
  354. FINSI ;
  355. FINSI ;
  356.  
  357. SI (> f 0.) ;
  358. *--- Valeurs propres relatives a lambda = i*omega = i*2*PI*f
  359. MESS i f (if/f) (-2*PI*if) (2*PI*f) Msg (((f**2) + (if**2))**0.5) (if/(((f**2) + (if**2))**0.5));
  360. SINON ;
  361. SI (EGA f 0.) ;
  362. MESS i f ' ------ ' ' ' (-2*PI*if) (2*PI*f) Msg (((f**2) + (if**2))**0.5) (if/(((f**2) + (if**2))**0.5)) ;
  363. FINSI ;
  364. FINSI ;
  365. Mod.Reel.i = Mod.Reel.i ET f ;
  366. Mod.Imag.i = Mod.Imag.i ET if ;
  367.  
  368. FIN Bloc5 ;
  369.  
  370. SI GRAPH ;
  371. MESS '-------------------------------------------------------------------------------' ;
  372. MESS ' ' ;
  373. FINSI ;
  374. FIN b_Vax ;
  375. *
  376. SAUTER 2 LIGNES ;
  377.  
  378. **********************************
  379. * Traces *
  380. **********************************
  381. * 4*2 couleurs
  382. numero = LECT Nbmod 'PAS' -1 1 1 'PAS' 1 Nbmod;
  383. coco4 = MOTS 'ROUG' 'ORAN' 'TURQ' 'BLEU' ;
  384. Nbcour = 2*Nbmod ;
  385. Tab1 = TABLE ;
  386. Tab2 = TABLE ;
  387. REPE Ident1 Nbcour ;
  388. * Tab1 . &Ident1 = 'NOLI MARQ S PLEIN TRIU' ;
  389. Tab1 . &Ident1 = 'NOLI MARQ S PLUS' ;
  390. Tab2 . &Ident1 = 'MARQ S PLUS' ;
  391. FIN Ident1 ;
  392.  
  393. *--- Traces relatifs a la frequence
  394. Evol1 = VIDE 'EVOLUTIO';
  395. Evol2 = VIDE 'EVOLUTIO';
  396. *--- Traces relatifs a OMEGA
  397. Evol3 = VIDE 'EVOLUTIO';
  398. Evol4 = VIDE 'EVOLUTIO';
  399. Evol5 = VIDE 'EVOLUTIO';
  400. * --- Boucle sur les modes
  401. REPE Trac1 Nbcour ;
  402. k = &Trac1 ;
  403. numk = EXTR numero k;
  404. si (k &lt;EG Nbmod); legek = CHAI 'Mode' ' ' numk;
  405. sinon; legek = MOT 'PAS DE LEGENDE';
  406. finsi;
  407. cocok = EXTR coco4 numk;
  408. SI (> (DIME Mod.Reel.k) 0) ;
  409. Evol1 = Evol1 ET (EVOL cocok MANU 'LEGE' legek
  410. 'V_{r}' LVR 'Re(f) (Hz)' (Mod.Reel.k)) ;
  411. Evol2 = Evol2 ET (EVOL cocok MANU 'LEGE' legek
  412. 'V_{r}' LVR 'Im(f) (/s)' (Mod.Imag.k)) ;
  413. Evol3 = Evol3 ET (EVOL cocok MANU 'LEGE' legek
  414. 'V_{r}' LVR 'Re(w)' (2.*PI*(Rw*Mod.Reel.k))) ;
  415. Evol4 = Evol4 ET (EVOL cocok MANU 'LEGE' legek
  416. 'V_{r}' LVR 'Im(w)' (2.*PI*Rw*(Mod.Imag.k))) ;
  417. Evol5 = Evol5 ET (EVOL cocok MANU 'LEGE' legek
  418. 'Re(w)' (2.*PI*Rw*(Mod.Reel.k)) 'Im(w)' (2.*PI*Rw*(Mod.Imag.k))) ;
  419. SINON;
  420. MESS 'Mode' k 'vide !';
  421. FINSI ;
  422. FIN Trac1 ;
  423.  
  424. SI GRAPH;
  425.  
  426. * avec legende
  427. DESS Evol1 'GRILL' 'GRIS' Tab2 'LEGE';
  428. * sans/avec zoom
  429. TITR 'Vitesse reduite - Frequences (partie reelle)' ;
  430. DESS Evol1 'GRILL' 'GRIS' Tab1 ;
  431. DESS Evol1 'XBOR' 0. 12. XGRA 1. 'YBOR' 0. 120. YGRA 10 'GRILL' 'GRIS' Tab1 ;
  432.  
  433. TITR 'Vitesse reduite - Frequences (partie imaginaire)' ;
  434. DESS Evol2 'GRILL' 'GRIS' Tab1 ;
  435. DESS Evol2 'XBOR' 0. 12. XGRA 1. 'YBOR' -0. 10. YGRA 1. 'GRILL' 'GRIS' Tab1 ;
  436.  
  437. TITR 'Frequences (partie reelle - partie imaginaire)' ;
  438. DESS Evol5 'GRILL' 'GRIS' Tab1 ;
  439. DESS Evol5 'XBOR' 0. 70. XGRA 10. 'YBOR' 0. 20. YGRA 2. 'GRILL' 'GRIS' Tab1 ;
  440. FINSI ;
  441.  
  442.  
  443.  
  444. *** TEST A PREVOIR ***
  445.  
  446.  
  447.  
  448. FIN ;
  449.  
  450.  
  451.  
  452.  
  453.  
  454.  
  455.  
  456.  
  457.  
  458.  
  459.  
  460.  
  461.  
  462.  
  463.  

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