* fichier : rotor5.dgibi
******************************************************************
*
*
* Etude d'une machine tournante
* dans le repere inertiel (ou fixe)
* Elements: POUT et volant d inertie
*
* Etude d'un arbre lors d'une rupture d'aube
*
* 1- Calcul des inerties et des forces de balourd a l aide
* d un modele 3D massif
* 2- Etude de la reponse a une rupture accidentelle d'aube
*
* Didier Combescure Juillet 2003
*
******************************************************************
*GRAPH = 'Y';
GRAPH = 'N';
*COMPLET = VRAI;
COMPLET = FAUX;
*
OPTI 'DIME' 3 'ELEM' SEG2;
*
* Caractéristique de l'arbre
*
L1 = 10.0;
Re1 = 0.2 ;
Re1S = 0.22 ;
E1 = 2.e+11 ;
Nu1 = 0.2 ;
Rho1 = 7800. ;
RhoFl = 1000.;
Visc1 = E1/(10000.);
Sarbre = pi*(Re1**2);
Inxarb = pi*(Re1**4)/2.;
Inyarb = pi*(Re1**4)/4.;
*
P0 = 0. 0. 0.;
P1 = (0.75*L1) 0. 0. ;
P2 = L1 0. 0. ;
*
Pb0 = 0. 0. 0.;
Pb2 = L1 0. 0.;
GeoRot = (D 75 P0 P1) et (D 25 P1 P2);
*
* Géométrie de la turbine
*
OPTI 'DIME' 3 ELEM CU20;
*
Re2 = 1.0;
Ep2 = 0.05;
vh2 = Ep2 0. 0.;
*
PM1= P1 plus (0. Re1 0.);
PM2 = P1 plus (0. 0. Re1);
PM3= P1 plus (0. ((-1.)*Re1) 0.);
PM4 = P1 plus (0. 0. ((-1.)*Re1));
PMrot = P1 plus (0. Re1 Re1);
*
PM1e= P1 plus (0. Re2 0.);
PM2e = P1 plus (0. 0. Re2);
PM3e= P1 plus (0. ((-1.)*Re2) 0.);
PM4e = P1 plus (0. 0. ((-1.)*Re2));
*
Lig2 = (CER3 4 PM1 PM2 PM3) ET (CER3 4 PM3 PM4 PM1);
Sur2 = SURF LIG2 'PLANE';
Vol2 = (VOLU Sur2 TRAN (0.5*VH2)) et (VOLU Sur2 TRAN ((-0.5)*VH2));
elim 0.00001 Vol2;
*
*
PM11e = d 2 PM1 PM1e;
Sur3a = PM11e ROTA 2 90. P0 P2;
Vol3a = (VOLU 1 Sur3a TRAN (0.5*VH2))
et (VOLU 1 Sur3a TRAN ((-0.5)*VH2));
DEPL Vol3a TOUR 25. P1 PMrot;
elim 0.001 Vol3a;
Vol3b = Vol3a TOUR 90. P0 P2;
Vol3c = Vol3b TOUR 90. P0 P2;
Vol3d = Vol3c TOUR 90. P0 P2;
VOL3 = Vol3a et Vol3b et Vol3c et Vol3d;
VOL3Bal = Vol3a et Vol3b et Vol3c;
*
* Modele de l'arbre
*
Mod1 = MODL GEOROT 'MECANIQUE' POUT;
Mat1 = MATE Mod1 'YOUN' E1 'NU ' Nu1 'RHO ' Rho1
'SECT' Sarbre 'INRY' Inyarb 'INRZ' Inyarb
'TORS' Inxarb 'OMEG' (1.0) 'VISQ' (Visc1);
*
MESH3D = POUT2MAS Mod1 Mat1 'MASSIF';
SI (EGA GRAPH 'Y');
OEIL = 100. 100. 100.;
ELIM 0.0001 MESH3D;
Titre 'Maillage poutre';
trac OEIL MESH3D CACH;
trac OEIL (MESH3D et GeoRot);
Titre 'Maillage poutre + volant';
trac OEIL (MESH3D et VOL3) CACH;
oeil = 100. 0. 0.;
Titre 'Volant';
trac oeil VOL3 cach;
trac oeil VOL3Bal cach;
FINSI;
*
* Modele de la turbine
*
Mod2 = MODL Vol3 'MECANIQUE' 'ELASTIQUE' ;
Mat2 = MATE Mod2 'YOUN' E1 'NU ' Nu1 'RHO ' Rho1;
*
Mod2B = MODL Vol3Bal 'MECANIQUE' 'ELASTIQUE' ;
Mat2B = MATE Mod2B 'YOUN' E1 'NU ' Nu1 'RHO ' Rho1;
*
*
Blo2 = (BLOQ 'UZ' P0) ET
(BLOQ 'UY' P0) ET
(BLOQ 'UX' P0);
Blo3 = (BLOQ 'UZ' P2) ET
(BLOQ 'UY' P2) ET
(BLOQ 'UX' P2);
Blo1 = BLOQ 'RX' Georot;
*
* Matrice de masse
*
MMass1 = MASS Mod1 Mat1 ;
MRigi1 = RIGI Mod1 Mat1 ;
*
* Condensation de la matrice de masse du volant d'inertie sur 1 point
*
MMass2 = MASS Mod2 Mat2 ;
MRigi2 = RELA CORI DEPL ROTA (P1 et Vol3);
*
SupK2 = SUPER 'RIGIDITE' 'NOMU' MRigi2 P1;
SupM2 = SUPER 'MASSE' SupK2 MMass2;
*
InxVol = EXTR (PSCAL (SupM2 *(MANU CHPO P1 1 RX 1.))
(MANU CHPO P1 1 RX 1.) (MOTS MX) (MOTS RX)) 'SCAL' P1;
InyVol = EXTR (PSCAL (SupM2 *(MANU CHPO P1 1 RY 1.))
(MANU CHPO P1 1 RY 1.) (MOTS MY) (MOTS RY)) 'SCAL' P1;
*
*
MMassT = MMass1 et SupM2;
*
Geovolan = MANU POI1 P1;
CG2 = MANU RIGI 'TYPE' 'AMORTISSEMENT' Geovolan
('MOTS' 'RX ' 'RY ' 'RZ ') 'QUEL'
(prog 0. 0. 0.
0. 0. (InyVol)
0. ((-1.)*InyVol) 0.);
*
* Matrice de couplage gyroscopique de l'arbre
*
CG1 = GYROSCOPIQUE Mod1 Mat1;
CGTOT = CG1 et CG2;
*
* Matrice d'amortissement
*
CDAM1 CDAH1 = 'AMOR' Mod1 Mat1 'COROTATIF';
*
*
* Calcul des forces de balourd (pour une vitesse de rotation de 1 tour/s)
*
* Définition du champs d'acceleration centrifuge
*
XYZ = (0.*(CHAN COMP (COOR 1 Vol3Bal) 'UX')) ET
(CHAN COMP (COOR 2 Vol3Bal) 'UY') ET
((CHAN COMP (COOR 3 Vol3Bal) 'UZ'));
*
XYZ0 = MANU CHPO Vol3Bal 3 'UX' (0.) 'UY' (0.)
'UZ' (0.);
*
* Calcul avec la matrice de raideur centrifuge
*
OmegRot = 1.0;
Vrota = OmegRot 0. 0.;
KCE2 ='KENT' 'CENTRIFUGE' Mod2b Mat2b Vrota;
FCENT2 = KCE2*(XYZ - XYZ0);
*
* Calcul avec la matrice de masse
*
FCENT = (OmegRot**2)*MMass2*(XYZ - XYZ0);
*
SI (EGA GRAPH 'Y');
oeil = 100. 0. 0.;
TRAC oeil Vol3Bal (vecteur FCENT FX FY FZ 0.01);
oeil = 0. 100. -100.;
TRAC oeil Vol3Bal (vecteur FCENT FX FY FZ 0.01);
FINSI;
*
* Condensation des forces centrifuges pour le calcul du balourd
*
Fbalourd = SUPER 'CHARGE' SupK2 FCENT;
FbalR = Fbalourd;
mess 'FbalR';
list FbalR;
Vecrot = MANU CHPO (extr FbalR mail) 3
RX 1. RY 0. RZ 0.;
FbalI = (FbalR '*'VecRot ('MOTS' FY FZ FX)
('MOTS' RZ RX RY) ('MOTS' FX FY FZ))
- (FbalR '*'VecRot ('MOTS' FZ FX FY)
('MOTS' RY RZ RX) ('MOTS' FX FY FZ));
FbalI = CHAN ATTRIBUT FbalI NATURE DISCRET;
**********************************
*
* Calcul des modes propres
*
**********************************
TBasR1 = VIBR 'PROCHE' (prog 10.) (lect 2)
(MRigi1 ET Blo1 ET Blo2 et Blo3)
MMassT 'TBAS' ;
*
TMod1 = TBasR1 . 'MODES' ;
NbMod1 = (DIME TMod1) - 2 ;
TInd1 = INDEX (TMod1 . 1) ;
*
SI (EGA GRAPH 'Y');
oeil = 100. 1. 0. ;
REPETER bbb1 NbMod1 ;
i = &bbb1 ;
Frq1 = TMod1 . i . 'FREQUENCE' ;
Def1 = TMod1 . i . 'DEFORMEE_MODALE' ;
TITRE Frq1 ;
trac oeil (defo georot def1);
FIN bbb1 ;
FINSI;
Frq1 = TMod1 . 1 . 'FREQUENCE' ;
*
*
MAmor1 = CDAM1;
*
*-Projection des matrices assemblées sur la base réelle
*
MMasP1 = PJBA TBasR1 MMassT ;
MRigP1 = PJBA TBasR1 (MRigi1 ET Blo1 ET Blo2 et Blo3) ;
MAmoP1 = PJBA TBasR1 MAmor1 ;
MAmHP1 = PJBA TBasR1 CDAH1;
Cg1P = PJBA TBasR1 CgTot ;
*
* Vitesse de rotation (10 Hz)
*
Omeg1 =10.*2.*pi;
*Omeg1 =0.*2.*pi;
*
* Calcul des modes complexes avec VIBC
*
BasCi1 = VIBC MMasP1 (MRigP1 et (Omeg1*MAmHP1))
(MAmoP1 et (Omeg1*Cg1P)) TBasR1 ;
*
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;
*
*
* Calcul du diagramme de Campbell
*
TAB1 = TABLE;
TAB1.'BASE_MODALE' = TBasR1;
TAB1.'MASS_PROJ' = MMasP1;
TAB1.'RIGI_PROJ' = MRigP1;
TAB1.'GYRO_PROJ' = 2.*pi*Cg1P;
TAB1.'AMOR_PROJ' = MAmoP1;
TAB1.'KROT_PROJ' = (2.*pi*MAmHP1);
TAB1.'AFFICHAGE' = VRAI;
SI COMPLET;
PRFREQ = (prog 0. pas 1. 50.);
SINON;
PRFREQ = prog 0. pas 1. 50.;
FINSI;
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;
*
EvFreOme = evol manu (prog 0. 50.) (prog 0. (50.));
EvFreZer = evol manu (prog 0. 50.) (prog 0. 0.);
*
SI (EGA GRAPH 'Y');
'DESS' (EVFRb1 et EVFRb2 et EvFreOme)
'TITR' 'Frequences' 'TITX' 'Freq rot (Tour/s)'
'TITY' 'Frequences (Hz)';
'DESS' (EVAMb1 et EVAMb2 et EvFreZer)
'TITR' 'Amortissement' 'TITX' 'Freq rot (Tour/s)'
'TITY' 'Amortissement (%)';
finsi;
*
SI COMPLET;
promeg = (prog 1. pas 0.1 50.)*2.*pi;
SINON;
promeg = (prog 1. pas 10. 50.)*2.*pi;
FINSI;
*
TAB1 = TABLE;
TAB1.'BASE_MODALE' = TBasR1 ;
TAB1.'MASS_PROJ' = MMasP1;
TAB1.'RIGI_PROJ' = MRigP1;
TAB1.'GYRO_PROJ' = Cg1P;
TAB1.'AMOR_PROJ' = (MAmoP1);
TAB1.'KROT_PROJ' = (MAmHP1);
TAB1.'AFFICHAGE' = VRAI;
TAB1.'FBALOURD' = Fbalourd;
TAB1.'SAUVDEFO' = FAUX;
TAB1.'REPONSE' = TABLE;
(TAB1.'REPONSE'). 1 = TABLE;
((TAB1.'REPONSE'). 1 ). 'POINT' = P1;
BALOURD TAB1 promeg;
*
EvUYbal = EVOL MANU (promeg/(2.*pi))
(extr (((TAB1.'REPONSE'). 1) . 'UY' ) ordo);
EvUZbal = EVOL MANU (promeg/(2.*pi))
(extr (((TAB1.'REPONSE'). 1) . 'UZ' ) ordo);
EvGabal = EVOL MANU (promeg/(2.*pi))
(extr (((TAB1.'REPONSE'). 1) . 'GAXE') ordo);
EvPabal = EVOL MANU (promeg/(2.*pi))
(extr (((TAB1.'REPONSE'). 1) . 'PAXE') ordo);
EvAngEl = EVOL MANU (promeg/(2.*pi))
(extr (((TAB1.'REPONSE'). 1) . 'ANGLEGAXE') ordo);
*
SI (EGA GRAPH 'Y');
'DESS' EvUYbal LOGY XBORD 0. 50.;
'FINSI';
**************************
*
* Calculs transitoires
*
**************************
*
* DYNAMIC
*
**************************
fr1 = 5.0;
dt = (1./(50.*fr1));
tdebut = (20.*dt);
SI COMPLET;
tduree = 20.*30.*dt;
SINON;
tduree = 3.*30.*dt;
FINSI;
*tduree = 0.25;
prtime2 = prog 0. pas dt (tduree + (100.*dt));
prtime = (prog 0.) et
(prtime2 + (prog (dime prtime2)*tdebut));
prcren = (prog 0.) et (prog (dime prtime2)*1.);
prcos = (prog 0. ) et (cos (fr1*prtime2*360.));
prsin = (prog 0. ) et (sin (fr1*prtime2*360.));
prtimef = prog 0. pas dt (tdebut+tduree);
EvReel = evol manu prtime
((prcos*prcren)*((2.*pi*fr1)**2));
EvImag = evol manu prtime
( (((-1.)*prsin*prcren)*((2.*pi*fr1)**2)));
*
ChaBal1 = (CHAR FbalR EvReel);
Chabal2 = (CHAR FbalI EvImag);
Chabal = Chabal1 et Chabal2;
*
LIS1 = prtimef;
NTT = ((dime prtimef) - 5);
TAB2 = TABLE ;
TAB2.'CHAR' = ChaBal ;
TAB2.'DEPL' = MANU CHPO GEOROT 3 UX 0. UY 0. UZ 0. ;
TAB2.'VITE' = MANU CHPO GEOROT 3 UX 0. UY 0. UZ 0. ;
TAB2.'RIGI' = (MRigi1 ET Blo1 ET Blo2 et Blo3)
et (0.*CDAH1*2.*pi*fr1);
TAB2.'MASS' = MMassT;
TAB2.'AMOR' = MAmor1 et (CgTot*2.*pi*fr1);
TAB2.'FREQ' = (1./dt);
TAB2.'INST' = LIS1 ;
TAB1 = DYNAMIC TAB2 ;
*
Time = PROG ;
PRUY1 = PROG ;
PRUZ1 = PROG ;
I = 0 ;
REPETER BB1 NTT ;
I = I + 1 ;
Time = Time ET ( PROG(TAB1 . I . 'TEMP' ));
DEP1 = TAB1 . I . 'DEPL' ;
PRUY1 = PRUY1 ET( PROG (EXTR DEP1 UY P1) );
PRUZ1 = PRUZ1 ET( PROG (EXTR DEP1 UZ P1) );
FIN BB1;
*
EvUY1 = evol manu Time PRUY1;
EvUZ1 = evol manu Time PRUZ1;
EvUYUZ1 = EVOL MANU PRUY1 PRUZ1;
*
SI (EGA GRAPH 'Y');
'DESS' EvUY1;
'DESS' EvUZ1;
'DESS' EvUYUZ1;
'FINSI';
*****************************
*
* DYNE
*
******************************
TRIG = TABLE 'RAIDEUR_ET_MASSE';
TRIG.'RAIDEUR' = MRigP1;
TRIG.'MASSE' = MMasP1;
*
TAMOR = TABLE 'AMORTISSEMENT';
TAMOR . 'AMORTISSEMENT' = MAmoP1;
*
*
CGP = (PJBA TBasR1 (CgTot*2.*pi*fr1));
*
TLIAI = TABLE 'LIAISON' ;
TTLA = TABLE 'LIAISON_A' ;
TLIAI.'LIAISON_A' = TTLA ;
*
* Couplage gyroscopique => mis au second membre
*
Nmode = (DIME (TBASR1.'MODES')) - 2;
I=0;
Iliai = 0;
REPETER BLOC1 Nmode;
I=I+1;
J=0;
REPETER BLOC2 Nmode;
J=J+1;
Iliai = Iliai + 1;
*
PI = ((TBASR1.'MODES'). I). 'POINT_REPERE';
PJ = ((TBASR1.'MODES'). J). 'POINT_REPERE';
*
Vunit = MANU CHPO PI 1 'ALFA' 1.;
CHPOM1 = CGP*Vunit;
CoefCG1 = EXTR CHPOM1 'FALF' PJ;
list CoefCG1;
TLIAI.'LIAISON_A'.Iliai=TABLE 'LIAISON_ELEMENTAIRE';
TLIAI.'LIAISON_A'.Iliai.'TYPE_LIAISON'= MOT 'COUPLAGE_VITESSE';
TLIAI.'LIAISON_A'.Iliai.'SUPPORT'=PJ;
TLIAI.'LIAISON_A'.Iliai.'ORIGINE'=PI;
TLIAI.'LIAISON_A'.Iliai.'COEFFICIENT'= (-1.)*CoefCG1;
*
FIN BLOC2;
FIN BLOC1;
*
TCHAR = TABLE 'CHARGEMENT' ;
CHAP1 = PJBA TBasR1 ChaBal;
TCHAR.'BASE_A' = CHAP1;
*
TSORT = TABLE 'SORTIE' ;
TSORV = TABLE 'VARIABLE' ;
TSORT.'VARIABLE' = TSORV ;
TSORV.'VITESSE' = VRAI;
TSORV.'DEPLACEMENT' = VRAI ;
TSORV.'ACCELERATION' = VRAI ;
*
DTEX = DT/2;
TRESU = DYNE DE_VOGELAERE TBasR1 TAMOR TCHAR TLIAI
(NTT*2) DTEX TSORT ;
************************************************
* TRACE DES DEPLACEMENTS EN FONCTION DU TEMPS
*************************************************
Ev2UY1 = EVOL VERT RECO TRESU TBasR1 'DEPL' P1 UY ;
Ev2UZ1 = EVOL BLEU RECO TRESU TBasR1 'DEPL' P1 UZ ;
Ev2UYUZ1 = EVOL ROUGE MANU (extr Ev2UY1 ordo) (extr Ev2UZ1 ordo);
*
SI (EGA GRAPH 'Y');
'DESS' (EvUY1 et Ev2UY1);
'DESS' (EvUZ1 et Ev2UZ1);
'DESS' (EvUYUZ1 et Ev2UYUZ1);
'FINSI';
*
UyMa= maxi (extr EvUY1 ordo) abs;
Uy2Ma= maxi (extr Ev2UY1 ordo) abs;
*
FR1Th = 7.496;
Uy2Th =1.820E-02;
UyTh =1.902E-02;
*
Errfrq = 0.01*FR1Th;
Erruy = 0.01*UyTh;
*
SI ((ABS(Frq1 - FR1Th) > Errfrq ) OU
(ABS(UyMa - UyTh) > Erruy ) OU
(ABS(Uy2Ma - Uy2Th) > Erruy ));
ERRE 5;
SINON;
ERRE 0;
FINSI;
*
FIN;
*
© Cast3M 2003 - All rights reserved.
Disclaimer