* fichier : rotor4.dgibi
******************************************************************
*
*
* Test de GYROSCOPIQUE, AMOR, CAMPBELL
* pour les elements de poutre
* Etude d'une machine tournante
* dans le repere inertiel (ou fixe)
* Elements: TUYAU et SECTION (fibre)
* POUT, TIMO en exemple
*
* 1-Influence des déformations de cisaillement
* 2-Etude de la stabilité avec un amortissement corotatif
*
* Modélisation d'un tube
*
* Didier Combescure Juillet 2003
*
******************************************************************
*GRAPH = 'Y';
GRAPH = 'N';
*
*L1 = 0.50 ;
L1 = 1.20;
Re1 = 0.1 ;
epa1 = 0.0012;
E1 = 2.e+11 ;
Nu1 = 0.2 ;
Rho1 = 1600. ;
Visc1 = 0.00001*E1;
*Visc1 = 0.0*E1;
*
* Description de la section
*
OPTI 'DIME' 2 'ELEM' QUA4;
P0S = 0. 0.;
P1S = (Re1 - (0.5*epa1)) 0.;
P2S = (Re1 + (0.5*epa1)) 0.;
*
P12S = D 1 P1S P2S;
SURS = P12S ROTA 50 360. P0S;
ELIM 0.0001 SURS;
*
OPTI 'DIME' 3 'ELEM' SEG2;
*
*
P0 = 0. 0. 0.;
P1 = L1 0. 0. ;
GeoRot = D 100 P0 P1;
*
* OMEG=1
* On donne une vitesse de rotation de 1 rad/s
*
* MODELE SECTION
*
kappa = 0.52;
*kappa = 1.0;
*kappB = 1.00;
*
MOSEC1 = MODE SURS MECANIQUE ELASTIQUE QUAS;
MASEC1 = MATE MOSEC1 'YOUNG' E1 'NU' Nu1 'RHO' Rho1
ALPY kappa ALPZ kappa 'VISQ' Visc1;
*
*
* Elément TIMO - FIBRE
*
Mod1 = MODL GEOROT 'MECANIQUE' ELASTIQUE SECTION TIMO;
Mat1 = MATR Mod1 'MODS' MOSEC1 'MATS' MASEC1
'OMEG' (2.*pi) VECT (0. 0. 1.) ;
*
*
* Elément TUYAU
*
*
Mod2 = MODL GEOROT 'MECANIQUE' TUYA;
Mat2 = MATE Mod2 'YOUNG' E1 'NU' Nu1 'RHO' Rho1
'RAYO' Re1 'EPAI' EPA1 'OMEG' (2.*pi) 'VISQ' Visc1;
*
* Elément POUT
*
*
Mod3 = MODL GEOROT 'MECANIQUE' POUT;
Mat3 = MATE Mod3 'YOUNG' E1 'NU' Nu1 'RHO' Rho1
'SECT' (pi * ((Re1 ** 2)-((Re1 - epa1)**2)) )
'INRY' (pi * ((Re1 ** 4)-((Re1 - epa1)**4))/4.)
'INRZ' (pi * ((Re1 ** 4)-((Re1 - epa1)**4))/4.)
'TORS' (pi * ((Re1 ** 4)-((Re1 - epa1)**4))/2.)
* 'SECY' (kappB*(pi * ((Re1 ** 2)-((Re1 - epa1)**2))) )
* 'SECZ' (kappB*(pi * ((Re1 ** 2)-((Re1 - epa1)**2))) )
'OMEG' (2.*pi) 'VISQ' Visc1;
*
*
* Timoshenko
*
Mod4 = MODL GEOROT 'MECANIQUE' TIMO;
Car4 = CARA Mod4 'SECT' (pi * ((Re1 ** 2)-((Re1 - epa1)**2)) )
'INRY' (pi * ((Re1 ** 4)-((Re1 - epa1)**4))/4.)
'INRZ' (pi * ((Re1 ** 4)-((Re1 - epa1)**4))/4.)
'TORS' (pi * ((Re1 ** 4)-((Re1 - epa1)**4))/2.)
'SECY' (kappa*(pi * ((Re1 ** 2)-((Re1 - epa1)**2))) )
'SECZ' (kappa*(pi * ((Re1 ** 2)-((Re1 - epa1)**2))) )
'OMEG' (2.*pi);
*
*
Mat4 = MATE Mod4 'YOUNG' E1 'NU' Nu1 'RHO' Rho1 'VISQ' Visc1;
Mat4 = Mat4 et Car4;
*
MESH3D = POUT2MAS Mod1 Mat1 'GAUSS';
SI (EGA GRAPH 'Y');
OEIL = 10. 10. 10.;
ELIM 0.0001 MESH3D;
Titre 'Maillage poutre avec section fibre';
trac OEIL MESH3D CACH;
trac OEIL (MESH3D et GeoRot);
FINSI;
*
*
* Cond. limites : 2 appuis
*
*
Blo2 = (BLOQ 'UZ' P0) ET
(BLOQ 'UY' P0) ET
(BLOQ 'UX' P0);
Blo3 = (BLOQ 'UZ' P1) ET
(BLOQ 'UY' P1) ET
(BLOQ 'UX' P1);
Blo1 = BLOQ 'RX' Georot;
*
MMass1 = MASS Mod1 (Mat1) ;
MRigi1 = RIGI Mod1 (Mat1) ;
*
MMass2 = MASS Mod2 (Mat2) ;
MRigi2 = RIGI Mod2 (Mat2) ;
*
MMass3 = MASS Mod3 (Mat3) ;
MRigi3 = RIGI Mod3 (Mat3) ;
*
MMass4 = MASS Mod4 (Mat4) ;
MRigi4 = RIGI Mod4 (Mat4) ;
*
* Matrice de couplage gyroscopique
*
CG1 = GYROSCOPIQUE Mod1 Mat1;
CG2 = GYROSCOPIQUE Mod2 Mat2;
CG3 = GYROSCOPIQUE Mod3 Mat3;
CG4 = GYROSCOPIQUE Mod4 Mat4;
*
CDAM1 CDAH1 = 'AMOR' Mod1 Mat1 'COROTATIF';
CDAM2 CDAH2 = 'AMOR' Mod2 Mat2 'COROTATIF';
CDAM3 CDAH3 = 'AMOR' Mod3 Mat3 'COROTATIF';
CDAM4 CDAH4 = 'AMOR' Mod4 Mat4 'COROTATIF';
*
TBasR1 = VIBR 'PROCHE' (prog 100.) (lect 2)
(MRigi1 ET Blo1 ET Blo2 et Blo3)
MMass1 'TBAS' ;
TBasR2 = VIBR 'PROCHE' (prog 100.) (lect 2)
(MRigi2 ET Blo1 ET Blo2 et Blo3)
MMass2 'TBAS' ;
TBasR3 = VIBR 'PROCHE' (prog 100.) (lect 2)
(MRigi3 ET Blo1 ET Blo2 et Blo3)
MMass3 'TBAS' ;
TBasR4 = VIBR 'PROCHE' (prog 100.) (lect 2)
(MRigi4 ET Blo1 ET Blo2 et Blo3)
MMass4 'TBAS' ;
*
TMod1 = TBasR1 . 'MODES' ;
TMod2 = TBasR2 . 'MODES' ;
TMod3 = TBasR3 . 'MODES' ;
TMod4 = TBasR4 . 'MODES' ;
NbMod3 = (DIME TMod3) - 2 ;
TInd3 = INDEX (TMod3 . 1) ;
*
SI (EGA GRAPH 'Y');
oeil = 100. 1. 0. ;
REPETER bbb1 NbMod3 ;
i = &bbb1 ;
Frq1 = TMod1 . i . 'FREQUENCE' ;
Def1 = TMod1 . i . 'DEFORMEE_MODALE' ;
Frq2 = TMod2 . i . 'FREQUENCE' ;
Def2 = TMod2 . i . 'DEFORMEE_MODALE' ;
Frq3 = TMod3 . i . 'FREQUENCE' ;
Def3 = TMod3 . i . 'DEFORMEE_MODALE' ;
Frq4 = TMod4 . i . 'FREQUENCE' ;
Def4 = TMod4 . i . 'DEFORMEE_MODALE' ;
TITRE Frq1 ;
trac oeil (defo georot def1);
TITRE Frq2 ;
trac oeil (defo georot def2);
TITRE Frq3 ;
trac oeil (defo georot def3);
TITRE Frq4 ;
trac oeil (defo georot def4);
FIN bbb1 ;
FINSI;
*
* On teste pour le modele a fibre
*
MAmor1 = CDAM1;
MAmor2 = CDAM2;
MAmor3 = CDAM3;
MAmor4 = CDAM4;
*
*-Projection des matrices assemblées sur la base réelle
*
MMasP1 = PJBA TBasR1 MMass1 ;
MRigP1 = PJBA TBasR1 MRigi1 ;
MAmoP1 = PJBA TBasR1 MAmor1 ;
MAmHP1 = PJBA TBasR1 CDAH1;
Cg1P = PJBA TBasR1 Cg1 ;
*
MMasP2 = PJBA TBasR2 MMass2 ;
MRigP2 = PJBA TBasR2 MRigi2 ;
MAmoP2 = PJBA TBasR2 MAmor2 ;
MAmHP2 = PJBA TBasR2 CDAH2;
Cg2P = PJBA TBasR2 Cg2 ;
*
MMasP3 = PJBA TBasR3 MMass3 ;
MRigP3 = PJBA TBasR3 MRigi3 ;
MAmoP3 = PJBA TBasR3 MAmor3 ;
MAmHP3 = PJBA TBasR3 CDAH3;
Cg3P = PJBA TBasR3 Cg3 ;
*
MMasP4 = PJBA TBasR4 MMass4 ;
MRigP4 = PJBA TBasR4 MRigi4 ;
MAmoP4 = PJBA TBasR4 MAmor4 ;
MAmHP4 = PJBA TBasR4 CDAH4;
Cg4P = PJBA TBasR4 Cg4 ;
*
* Vitesse de rotation (0 Hz)
*
Omeg1 =0.;
*
* Calcul des modes complexes avec VIBC - Bernoulli
*
BasCi1 = VIBC MMasP4 (MRigP4 et (Omeg1*MAmHP4))
(MAmoP4 et (Omeg1*Cg4P)) TBasR4 ;
*
ModBasC = BasCi1.'MODES';
NbMod1 = (DIME ModBasC) - 2;
*
FRQR = TABLE;
FRQI = TABLE;
FRQ = TABLE;
AMOR= TABLE;
*
mess 'Frequences et amortissement au repos';
*
REPETER bbb2 NbMod1 ;
i2 = &bbb2 ;
FrqR1 = BasCi1.'MODES'.i2.'FREQUENCE_REELLE' ;
FrqI1 = BasCi1.'MODES'.i2.'FREQUENCE_IMAGINAIRE' ;
DefR1 = BasCi1.'MODES'.i2.'DEFORMEE_MODALE_REELLE';
DefI1 = BasCi1.'MODES'.i2.'DEFORMEE_MODALE_IMAGINAIRE';
FRQR.i2 = FrqR1;
FRQI.i2 = FrqI1;
FRQ.i2 = ((FrqR1 ** 2) + (FrqI1 ** 2)) ** 0.5 ;
AMOR.i2 = FrqI1 /(FRQ.i2 );
mess 'Valeur propre réelle';
list (FrqR.i2);
mess 'Valeur propre imaginaire';
list (FrqI.i2);
mess 'Fréquence';
list (FRQ.i2);
mess 'Amortissement';
list (AMOR.i2);
FIN bbb2;
*
*
* Vitesse de rotation (1000 Hz)
*
Omeg1 =1000.;
*
* Calcul des modes complexes avec VIBC - Timoshenko
*
BasCi1 = VIBC MMasP4 (MRigP4 et (Omeg1*MAmHP4))
(MAmoP4 et (Omeg1*Cg4P)) TBasR4;
*
ModBasC = BasCi1.'MODES';
NbMod1 = (DIME ModBasC) - 2;
*
FRQR = TABLE;
FRQI = TABLE;
FRQ = TABLE;
AMOR= TABLE;
*
mess 'Frequences et amortissement a 1000 Hz';
*
REPETER bbb2 NbMod1 ;
i2 = &bbb2 ;
FrqR1 = BasCi1.'MODES'.i2.'FREQUENCE_REELLE' ;
FrqI1 = BasCi1.'MODES'.i2.'FREQUENCE_IMAGINAIRE' ;
DefR1 = BasCi1.'MODES'.i2.'DEFORMEE_MODALE_REELLE';
DefI1 = BasCi1.'MODES'.i2.'DEFORMEE_MODALE_IMAGINAIRE';
FRQR.i2 = FrqR1;
FRQI.i2 = FrqI1;
FRQ.i2 = ((FrqR1 ** 2) + (FrqI1 ** 2)) ** 0.5 ;
AMOR.i2 = FrqI1 /(FRQ.i2);
mess 'Valeur propre réelle';
list (FrqR.i2);
mess 'Valeur propre imaginaire';
list (FrqI.i2);
mess 'Fréquence';
list (FRQ.i2);
mess 'Amortissement';
list (AMOR.i2);
FIN bbb2;
*
* Calcul du diagramme de Campbell - Bernoulli
*
TAB1 = TABLE;
TAB1.'BASE_MODALE' = TBasR3;
TAB1.'MASS_PROJ' = MMasP3;
TAB1.'RIGI_PROJ' = MRigP3;
TAB1.'GYRO_PROJ' = Cg3P;
TAB1.'AMOR_PROJ' = MAmoP3;
TAB1.'KROT_PROJ' = MAmHP3;
TAB1.'AFFICHAGE' = VRAI;
*PRFREQ = (prog 0. pas 2.5 100.)*100/(2.*pi);
PRFREQ = (prog 0. pas 25 100.)*100/(2.*pi);
CAMPBELL TAB1 PRFREQ;
*
EVFRb1 = (TAB1. 1). 'FREQUENCE_MODULE' ;
EVFRb2 = (TAB1. 4). 'FREQUENCE_MODULE' ;
EVAMb1 = (TAB1. 1). 'AMORTISSEMENT' ;
EVAMb2 = (TAB1. 4). 'AMORTISSEMENT' ;
*
SI (EGA GRAPH 'Y');
dess (EVFRb1 et EVFRb2 );
dess (EVAMb1 et EVAMb2);
FINSI;
*
* Calcul du diagramme de Campbell - Timoshenko
*
TAB1 = TABLE;
TAB1.'BASE_MODALE' = TBasR4 ;
TAB1.'MASS_PROJ' = MMasP4;
TAB1.'RIGI_PROJ' = MRigP4;
TAB1.'GYRO_PROJ' = Cg4P;
TAB1.'AMOR_PROJ' = MAmoP4;
TAB1.'KROT_PROJ' = MAmHP4;
TAB1.'AFFICHAGE' = VRAI;
*PRFREQ = (prog 0. pas 2.5 100.)*100/(2.*pi);
PRFREQ = (prog 0. pas 25. 100.)*100/(2.*pi);
CAMPBELL TAB1 PRFREQ;
*
EVFRT1 = ((TAB1. 1). 'FREQUENCE_MODULE') COUL rouge ;
EVFRT2 = ((TAB1. 4). 'FREQUENCE_MODULE') COUL rouge ;
EVAMT1 = ((TAB1. 1). 'AMORTISSEMENT') COUL rouge ;
EVAMT2 = ((TAB1. 4). 'AMORTISSEMENT') COUL rouge ;
*
EvFreOme = evol manu (prog 0. 1600.) (prog 0. (1600.));
EvFreZer = evol manu (prog 0. 1600.) (prog 0. 0.);
*
SI (EGA GRAPH 'Y');
'DESS' (EVFRT1 et EVFRT2 et EvFreOme)
'TITR' 'Frequences' 'TITX' 'Freq rot (Tour/s)'
'TITY' 'Frequences (Hz)';
'DESS' (EVFRT1 et EVFRT2)
'TITR' 'Frequences' 'TITX' 'Freq rot (Tour/s)'
'TITY' 'Frequences (Hz)';
'DESS' (EVAMT1 et EVAMT2 et EvFreZer)
'TITR' 'Amortissement' 'TITX' 'Freq rot (Tour/s)'
'TITY' 'Amortissement (%)';
*
dess (EVFRT1 et EVFRT2 et EVFRB1 et EVFRB2 et EvFreOme);
dess (EVAMT1 et EVAMT2 et EVAMB1 et EVAMB2);
*
finsi;
*
mess 'Frequences à 0 Hz et 1000 Hz';
mess 'Poutre de Bernoulli';
FR1b = extr (IPOL (prog 0.01 1000.) (extr EVFRB1 absc)
(extr EVFRB2 ordo)) 1;
AM1b = extr (IPOL (prog 0.01 1000.) (extr EVAMB1 absc)
(extr EVAMB1 ordo)) 1;
FR2b = extr (IPOL (prog 0.01 1000.) (extr EVFRB1 absc)
(extr EVFRB2 ordo)) 2;
AM2b = extr (IPOL (prog 0.01 1000.) (extr EVAMB1 absc)
(extr EVAMB1 ordo)) 2;
mess 'Poutre de Timoshenko';
FR1T = extr (IPOL (prog 0.01 1000.) (extr EVFRT1 absc)
(extr EVFRT2 ordo)) 1;
AM1T = extr (IPOL (prog 0.01 1000.) (extr EVAMT1 absc)
(extr EVAMT1 ordo)) 1;
FR2T = extr (IPOL (prog 0.01 1000.) (extr EVFRT1 absc)
(extr EVFRT2 ordo)) 2;
AM2T = extr (IPOL (prog 0.01 1000.) (extr EVAMT1 absc)
(extr EVAMT1 ordo)) 2;
*
* On norme
*
EVFRb1N = EVOL MANU (extr EVFRb1 absc)
((extr EVFRb1 ordo)/FR1b );
EVFRb2N = EVOL MANU (extr EVFRb1 absc)
((extr EVFRb2 ordo)/FR1b );
EVFRT1N = EVOL MANU (extr EVFRT1 absc)
((extr EVFRT1 ordo)/FR1T );
EVFRT2N = EVOL MANU (extr EVFRT1 absc)
((extr EVFRT2 ordo)/FR1T );
*
SI (EGA GRAPH 'Y');
TA= table;
ta.1 = 'TIRR';
ta.2 = 'TIRR';
titre 'Diagramme de Campbell
(Frequence = f(Fréquence rotation en Hz)';
DESS (EVFRb1 et EVFRb2 et EVFRT1 et EVFRT2) ta ybord 0. 4000.;
DESS (EVFRb1 et EVFRb2 et EVFRT1 et EVFRT2) ta ;
titre 'Diagramme de Campbell
(F/F0 = f(Fréquence rotation en Hz)';
DESS (EVFRb1N et EVFRb2N et EVFRT1N et EVFRT2N) ta ybord 0. 2.;
DESS (EVFRb1N et EVFRb2N et EVFRT1N et EVFRT2N) ta;
FINSI;
*
*
FR1Th = 787.3;
FR2Th = 812.7;
AM1Th = 2.48E-02;
AM2Th = -0.559E-02;
Errfrq = 0.01*FR1Th;
Erramo = 0.01*AM1Th;
*
SI ((ABS(FR1T - FR1Th) > Errfrq ) OU
(ABS(FR2T - FR2Th) > Errfrq ) OU
(ABS(AM1T - AM1Th) > Erramo ) OU
(ABS(AM2T - AM2Th) > Erramo ));
ERRE 5;
SINON;
ERRE 0;
FINSI;
*
FIN;
*
© Cast3M 2003 - All rights reserved.
Disclaimer