Levent MALGACA*, Hira KARAGÜLLE*
*Dokuz Eylül Üniversitesi, Mühendislik Fakültesi, Makine Mühendisliği Bölümü 35100, İZMİR levent.malgaca@deu.edu.tr, hira.karagulle@deu.edu.tr
ÖZET
Bu çalışmada 3 serbestlik dereceli kütle yay sisteminin aktif titreşim kontrolü üç farklı yöntemle ele alınmıştır. İlk yöntemde, hareket denklemleri elde edildikten sonra analitik çözüm Laplace yöntemi ile bulunur. İkinci yöntemde, ANSYS programında sonlu elemanlar modeli ve aktif titreşim kontrolü gerçekleştirilir. Son yöntemde ise, çözüm entegre yaklaşım ile yapılır; ANSYS de kurulan modelin sonlu eleman matrisleri elde edildikten sonra kontrol kısmı MATLAB/Simulink programında gerçekleştirilir. Sonlu eleman matrisleri ANSYS programında bir çıktı dosyasına yazdırılır. Harwell- Boeing biçimindeki bu matrisler MATLAB programı ile kütle ve rijitlik matrislerine dönüştürülür.
Sistemin dinamik denklemleri durum değişkenleri biçiminde düzenlenir. Aktif kontrolün etkisini göstermek amacıyla sönüm ihmal edilir. Performans kriteri olarak uç noktadaki son kütlenin yer değiştirmesinin dinamik cevabı PID kontrolün değişik katsayıları için değerlendirilir. Zeminden adım girdiyle zorlanan sistem için, üç yöntemin sonuçları birbirleriyle karşılaştırılır. Bu çalışmada dikkate alınan üç yöntemin yaklaşık olarak aynı sonuçları verdiği gözlemlenir. ANSYS gibi sonlu elemanlar programlarına kontrol algoritmaları entegre edilerek aktif titreşim kontrolü gerçekleştirilir. Entegre analiz ile elde edilen sonuçlar, karmaşık sistemlerin sonlu eleman modelleriyle kontrol algoritmalarının bütünleştirilebileceğini ortaya çıkarır.
Anahtar Kelimeler: Aktif titreşim kontrolü, sonlu eleman, ANSYS, MATLAB/Simulink
INTEGRATION OF ACTIVE VIBRATION CONTROL WITH FINITE ELEMENT MODELS OF MECHANICAL SYTEMS
ABSTRACT
In this study, active vibration control of a 3-DOF spring-mass system is considered with three different methods. In the first method, analytical solution is found by the Laplace method after the equations of motion for the system are obtained. In the second method, finite element modeling and active vibration control are performed in ANSYS. In the last method, solution is done by the integrated approach; the control part is realized in MATLAB/Simulink after the finite element matrices of the system modeled by ANSYS are obtained. Finite element matrices are written to an output file in ANSYS. These matrices in Harwell-Boeing format are converted to mass and stiffness matrices by a MATLAB program. Dynamic equations of the system are arranged by the state space variables format.
Damping is ignored to show the effect of the active control. The dynamic response of the displacement of the last mass at the tip is evaluated as performance criteria for various gains of PID control. The results of the three methods are compared to each other for the system excited by a step input from the base. It is observed that the three methods considered in this study give approximately the same results. Active vibration control is realized by integrating control algorithms to the finite element programs like ANSYS. The results obtained by the integrated analysis reveal that finite element models and control algorithms can be integrated.
Keywords: Active vibration control, finite element, ANSYS, MATLAB/Simulink
1. GİRİŞ
Aktif titreşim kontrolü makina mühendisliği, uçak mühendisliği ve inşaat mühendisliği gibi değişik mühendislik sahalarında farklı yapılara uygulanabilmektedir. Uçak mühendisliğinde bir uçağın kanat veya kuyruğunun, inşaat mühendisliğinde bir binanın veya köprünün, makina mühendisliğinde ise robot kollarının veya esnek kirişlerin ve plakaların aktif titreşim kontrolü, mühendislerin ilgisini çeken araştırma konularıdır [1-6]. Esnek yapılar piezoelektrik aktüatörler ile, daha ağır yapılar ise servo-hidrolik elemanlar ile kumanda edilerek aktif kontrol başarılabilmektedir.
Kuvvet ve ivme gibi farklı geri beslemeler kullanılarak da aktif titreşim izolasyonu yapılabilmektedir [7].
O’Connor ve Lang [8] tek uzuvlu esnek bir robot kolunu topaklanmış kütleler ve yaylardan oluşan çok serbestlik dereceli bir sistem olarak modellemiştir. Aktif kontrol dalga yayılımı yöntemi ile başarılmıştır. Böylece esnek robot kollarının pozisyon kontrolü minimum artık titreşim ve daha hızlı zaman cevabı ile sağlanır.
Yazarlar önceki çalışmasında mekanik sistemlerin sonlu eleman modelleri ile kontrol algoritmalarını bütünleştirmişlerdir [9].
ANSYS programını kullanarak sonlu eleman analizine kontrol algoritmalarını dahil ederek aktif titreşim kontolünü ele almışlardır. Önce iki serbestlik dereceli mekanik bir sistemin aktif kontrolü analitik yöntem ile ele alınmış ve bulunan sonuçlar ANSYS çözümleri ile karşılaştırılmıştır. Doğrulama yapıldıktan sonra piezoelektrik aktüatör içeren akıllı bir kirişin aktif titreşim kontrolü ANSYS programı ile gerçekleştirilmiştir.
Bu çalışmada mekanik sistemlerin sonlu elemanlar modelleri ile aktif titreşim kontrolünün entegrasyonu ele alınmıştır.
Mekanik sistem olarak 3 serbestlik dereceli bir kütle yay titreşim sistemi dikkate alınmıştır.
Analitik, ANSYS ve entegre yaklaşım çözümleri üç farklı yöntem ile gerçekleştirilerek bulunanan sonuçlar sunulmuştur. Entegre yaklaşımda sonlu elemanlar matrisleri ANSYS den elde edilerek, kontrol kısmı MATLAB/Simulink’ de yerine getirilmiştir. Sonuçların birbiri ile uyumlu olduğu gözlemlenmiştir.
2. AKTİF TİTREŞİM KONTROLÜ
Ele alınan 3 serbestlik dereceli kütle-yay sistemi Şekil 1’ de gösterilmiştir. Burada m1, m2, m3 ve k1, k2, k3 kütle ve yay için parametrelerdir. Sistemin girdisi z(t), çıktıları ise her bir kütlenin yerdeğiştirmeleri x1, x2 ve x3 tür. Sistem zeminden birim adım girdi ile uyarılmaktadır. Aktif titreşim kontrolü uygulanarak uç noktadaki son kütlenin zamana göre yer değiştirmesi x3(t) sistem cevabı olarak değerlendirilmektedir. Sonuçlar üç farklı yöntem için de aynı oransal (Kp), integral (KI) ve türevsel (KD) kontrol kazançları ile elde edilir.
Şekil 1. 3 serbestlik dereceli kütle-yay sistemi (m1=m2=m3=1 kg, k1=k2=k3=300 N/m) 2.1 Analitik çözüm
Aktif titreşim kontrolünün analizi için ele alınan sistemin blok diyagramı Şekil 2’ de gösterilmiştir. Çok serbestlik dereceli sönümsüz bir titreşim sistemi için hareket denklemi denklem (1) deki gibi ifade edilebilir.
Burada M ve K sırasıyla, kütle ve rijitlik matrisleri, x değişkeni genel koordinat, f ise genel kuvvettir.
K
Şekil 2. 3 serbestlik dereceli kütle-yay sisteminin blok diyagramı
[ ] [ ]
M •x•+ Kx={ }
f (1) Ele alınan sistem için kinetik ve potansiyel enerjiler yazıldıktan sonra Lagrange denklemi uygulanarak hareket denklemleri matris biçiminde elde edilir.) t ( 0 z 0 k x x x k k 0
k k k k
0 k k k
x x x
m 0 0
0 m 0
0 0
m 1
3 2 1 3 3
3 3 2 2
2 2 1
3 2 1
3 2 1
⎪⎭
⎪⎬
⎫
⎪⎩
⎪⎨
=⎧
⎪⎭
⎪⎬
⎫
⎪⎩
⎪⎨
⎧
⎥⎥
⎦
⎤
⎢⎢
⎣
⎡
−+ −
−+ − +
⎪⎪
⎪
⎭
⎪⎪⎪
⎬
⎫
⎪⎪
⎪
⎩
⎪⎪⎪
⎨
⎧
⎥⎥
⎦
⎤
⎢⎢
⎣
⎡
•
•
•
•
•
•
(2)
Eğer x(t)=Xiest, z(t)=0 ise, sistemin sönümsüz doğal frekansları aşağıdaki gibi bulunabilir.
⎪⎭
⎪⎬
⎫
⎪⎩
⎪⎨
=⎧
⎪⎭
⎪⎬
⎫
⎪⎩
⎪⎨
⎧
⎥⎥
⎥
⎦
⎤
⎢⎢
⎢
⎣
⎡
+
−
− + +
−
− +
+
00 0 X X X k s m k 0
k k k s m k
0 k
) k k ( s m
3 2 1 2 3 3 3
3 3 2 2 2 2
2 2
2 1
1 (3)
Sayısal değerler yerine konarak A sistem matrisi aşağıdaki gibi elde edilir.
⎪⎭
⎪⎬
⎫
⎪⎩
⎪⎨
=⎧
⎪⎭
⎪⎬
⎫
⎪⎩
⎪⎨
⎧
⎥⎥
⎥
⎦
⎤
⎢⎢
⎢
⎣
⎡
+
−+ −
−
− +
00 0 X XX 300 s 300 0
300 600 s 300
0 300 600 s
3 2 1
A 2 2 2
4 4 4 4
4 3
4 4 4 4
4 2
1
(4)
A matrisinin determinantını sıfır yapan değerler sistemin özdeğerleridir.
s6+1500s4+540000s2+27000000=0 (5) ω1 = 7.7084 (rad/s) → f1 = 1.2268 (Hz) ω2 = 21.5983 (rad/s) → f2 = 3.4375 (Hz) ω3 = 31.2105 (rad/s) → f3 = 4.9673 (Hz)
Eğer girdi z(t)= est ve cevap x(t)=H(s) est ise, bu dönüşüm ile sistemdeki son kütlenin girdisi ve çıktısı arasındaki ilişkiyi gösteren açık devre transfer fonksiyonu Laplace yöntemi ile aşağıdaki denklemden elde edilebilir.
⎪⎭
⎪⎬
⎫
⎪⎩
⎪⎨
⎧
⎥⎥
⎥⎥
⎦
⎤
⎢⎢
⎢⎢
⎣
⎡
+ + +
+ +
+ +
+ +
+
⎪⎭=
⎪⎬
⎫
⎪⎩
⎪⎨
⎧
00 ) s ( 1Z k
270000 s2 4 1200 s ) 2 600 s ( 300 90000
) 2 600 s ( 300 180000 s2 4 900 s ) 2 300 s ( 300
90000 )
2 300 s ( 300 90000 s2 4 900 s
) s ( D
1 ) s 3( X2(s) X
) s 1( X
(6)
) s ( D 27000000 )
s ( Z
) s ( 3 ) X s (
HOL = = (7)
Benzer şekilde her bir kütle ile girdi arasındaki açık devre transfer fonksiyonları yazıldıktan sonra sistemin kapalı devre transfer fonksiyonu bulunabilir. Transfer fonksiyonlarının (H11, H12, H13) alt indislerinde; birinci indis girdiyi, ikinci indis ise her bir kütlenin genel koordinatını belirtmektedir.
X1(s) = H11(s).F(s) X2(s) = H12(s).F(s) X3(s) = H13(s).F(s)
F(s)=k1Z(s), Z(s)=1/s
) s ( D
90000 s
4 900 ) s s ( H
2
11 = + + (8)
) s ( D
90000 s
) 300 s ( H
2
12 = + (9)
) s ( D 90000 )
s (
H13 = (10)
Böylece sistemin blok diyagramından faydalanılarak kapalı devre transfer fonksiyonu hesaplanabilir. Burada referans değeri 1 olarak alınır.
U(s)={Xr(s)-K[X1(s)+X2(s)+X3(s)]} (11) X1(s) = uH11(s)
X2(s) = uH12(s) X3 (s)= uH13(s)
[H (s) H (s) H (s)]
) s ( H
) s ( )X s ( KG ) s ( G ) s ( ) x s ( H
) s ( X
13 12 11 13 1 3 1 r 13
3 = − + + (12)
Referans değer ile uç noktadaki kütlenin yer değiştirmesi arasındaki kapalı devre transfer fonksiyonu aşağıdaki gibi ifade edilir.
[H (s) H (s) H (s)]
) s ( KG 1
) s ( G ) s ( H )
s ( X
) s ( ) X s ( H
13 12 11 1
1 13 r
CL = 3 = + + + (13)
Sayısal değerler yerine konarak kapalı devre transfer fonksiyonu aşağıdaki gibi bulunur.
6 i p 6
6 i 3
p 4 3
4 i d 4 5
p 6
d 7
i p 2 d 6 CL
K 10 x 27 s ) 10 x 27 K 106 x 27 (
s ) 10 x 27 K 104 x 12 ( s ) K 10 x 12 10 x 540 (
s ) K 100 K 10 x 12 ( s ) K 100 1500 ( s K 100 s
) K s K s K ( 10 x ) 27
s ( H
2 +
+ +
+ +
+ +
+ +
+ + +
+
= +
(14)
Uç noktadaki son kütlenin zaman cevabı x3(t) denklem (14)’ de ters Laplace dönüşümü yapılarak bulunabilir. Analitik çözüm ile elde edilen kontrol edilmemiş ve kontrol edilmiş zaman cevapları Şekil 3(a) da gösterilmiştir.
2.2 ANSYS ile çözüm
Sistemin sonlu elemanlar modeli Şekil 1’de gösterilen düğüm numaralarına göre kurulur.
Sonlu eleman modelini kurmak için MASS21 ve COMBIN14 elemanlarından faydalanılır.
MASS21 topaklanmış kütle, COMBIN14 ise sönüm ve yay bilgisini içerir. Sistemde sönüm olmadığı için sönüm değeri sıfır alınarak sadece yaylar için değerler tanımlanır. Modal analiz yapılarak sistemin sönümsüz doğal frekansları bulunur. Geçici rejim (transient) analizi ile kapalı kontrol yapılır. Aktif kontrol aşağıdaki makro program ile yine ANSYS programı içinde gerçekleştirilir [9].
sum=0 errp=0
*do,t,2*dt,tson,dt
*get,x1,node,1,u,x
*get,x2,node,2,u,x
*get,x3,node,3,u,x err=1-(x1+x2+x3)/3 sum=sum+err*dt diff=(err-errp)/dt
ucon=kp*err+ki*sum+kd*diff d,4,ux,ucon
errp=err time,t solve
*enddo
Burada dt zaman adımıdır ve en küçük sönümsüz doğal frekansın (f1) 20 de biri alınabilir. Fakat türevsel kontrolün etkinliğini arttırabilmek amacıyla zaman adımı daha küçültülerek dt=1/f1/60=0.00335 s olarak alınmıştır. Düzenli rejim zamanı (tson) 5 s olarak alınmıştır. Kapalı kontrol “*do-*enddo”
döngüsünde gerçekleştirilir. Döngüden hemen önce zemini temsil eden 4 nolu düğüm noktasına birim adım girdi uygulanır ve ilk zaman adımı için çözüm yapılır. Daha sonra ikinci adımda, kütlelerin yerdeğişimi ilk adımdan bilindiği için kontrol döngüsü içinde hata sinyali hesaplanır. Seçilen kontrol parametrelerinin değerleri için kontrol yer değiştirmesi hesaplanarak aktüatörün varolduğu kabul edilen zemindeki düğüm noktasına uygulanır. Böylece çözüm düzenli
rejim zamanına kadar devam ettirilir. Daha iyi kontrol için geri beslemede yer değiştirmelerin aritmetik ortalaması dikkate alınmıştır.
ANSYS ile bulunan sonuçlar Şekil 3(b) de gösterilmiştir.
Şekil 3. Uç noktadaki son kütlenin yer değiştir mesinin kontrol edilmiş (Kp=4, KI=1, KD=0.1) ve kontrol edilmemiş dinamik cevapları (a) Analitik çözüm, (b) ANSYS ile çözüm,
(c) Entegre yaklaşım ile çözüm.
2.3 Entegre yaklaşım ile çözüm
Entegre analiz ANSYS, VisualBASIC ve MATLAB/Simulink programları kullanılarak gerçekleştirilir. Bölüm 2.2 de kurulan sonlu eleman modelinden faydalanılır. Sistemin kütle ve rijitlik matrisleri elde edildikten sonra aktif kontrol MATLAB da geliştirilen bir program yardımıyla yapılır. Entegre analiz için akış şeması Şekil 4’ de gösterilmiştir.
Şekil 4. Entegre yaklaşım ile aktif titreşim kontrolünün akış şeması
Akış şemasında da gösterildiği gibi sonlu eleman modeli kurulur ve daha sonra modal analiz yapılır. Böylece sonlu eleman sistem matrisleri “Harwell-Boeing” adı verilen kapalı biçimde bir çıktı dosyasına yazdırılır. Kütle ve rijitlik matrislerini içeren çıktı dosyaları “full”
uzantılıdır. VisualBASIC’ de geliştirilen bir program (kmmat.exe) ile “full” uzantılı dosyaların içindeki sayısal değerler okunarak
“m” uzantılı dosyalara yazdırılır.
Şekil 5. MATLAB/Simulink ile sistemin açık devre blok diyagramı
Şekil 6. Tüm durum değişkenlerinin açık devre zaman cevapları
Şekil 7. MATLAB/Simulink ile sistemin kapalı devre blok diyagramı
Şekil 8. Tüm durum değişkenlerinin kapalı devre zaman cevapları
MATLAB’ da geliştirilen bir program (cansys.m) ile kapalı biçimde bulunan kütle ve rijitlik matrisleri, bilinen biçimdeki global sistem matrislerine dönüştürülür. Bu aşamadan sonra sistemin aktif titreşim kontrolüne geçilir.
Sistemin dinamik denklemleri yine MATLAB’
da geliştirilen bir program (pautoss.m) ile durum değişkenleri biçimine getirilir. Son
olarak aktif titreşim kontrolü MATLAB/Simulink’ de uygulanr.
Açık devre blok diyagramı ve tüm durum değişkenlerinin zaman cevapları Şekil 5 ve 6’
da gösterilmiştir. Aktif titreşim kontrolünün gerçekleştirildiği kapalı devre blok diyagramı ile durum değişkenlerinin zaman cevapları Şekil 7 ve 8’ de gösterilmiştir. Diğer yöntemler ile karşılaştırabilmek için son kütlenin kontrol edilmemiş ve kontrol edilmiş zaman cevapları tek bir grafik üzerinde Şekil 3 (c)’ de gösterilmiştir.
3. SONUÇLAR
Basit geometrilere sahip yapıların sonlu eleman modelleri, bilinen formülasyonları kullanarak kurulabilmektedir. Daha karmaşık ve içinde farklı elemanları barındıran yapıların sonlu elemanlar modelleri bilgisayar destekli mühendislik programları yardımıyla kurulabilir. Bu çalışmada bahsedilen entegre yaklaşım ile sonlu eleman matrisleri elde edilebilir. Daha sonra aktif kontrol başka bir programda gerçekleştirilebilir. Burada entegre yaklaşım 3 serbestlik dereceli kütle-yay sistemine uygulanmıştır. Entegre yaklaşım ile elde edilen sonuçların analitik ve ANSYS sonuçlarıyla yaklaşık olarak aynı olduğu gözlemlenmiştir. Bu yaklaşım, karmaşık sistemlerin sonlu elemanlar modellerinin daha hızlı bir şekilde kurulabilmesine ve farklı kontrol algoritmaları ile bütünleştirilmesine imkan tanımaktadır.
4. KAYNAKLAR
1. Ülker F D, Yaman Y, Nalbantoğlu V, Kırcalı Ö F, Akıllı bir plağın burulma titreşimlerinin aktif kontrolü, 12.
Ulusal Mak Teo Semp, Bil Kitabı, 1 177- 182,9-11 Haziran 2005, Kayseri, Türkiye.
2. Seto K, Matsumoto Y, Active vibration control of multiple buildings connected with active control bridges in response to large earthquakes, Proceedings of the American Control Conference, 1007-1011, San Diego, California June 1999.
3. Reaves M C, Horta G H, Piezoelectric actuator modeling using MSC/NASTRAN and MATLAB, Langley Research Center, Hampton, Virginia, NASA/TM (2003) 212651.
4. Xu S X, Koko T S, Finite element analysis and design of actively controlled piezoelectric smart structures, Finite elements in analysis and design, 40 (2004) 241-262.
5. Xianmin Z, Changjian S, Erdman A G, Active vibration controller design and comparison study of flexible linkage mechanism systems, Mechanisms and Machine Theory 37 (2002) 985-997.
6. Sun D, Mills J K, Shan J, Tsoa S K, A PZT actuator control of a single-link flexible manipulator based on linear velocity feedback and actuator placement, Mechatronics, 14 (2004) 381-401.
7. Preumont A, François A, Bossens F, Hanieh A A, Force feedback versus acceleration feedback in active vibration isolation, Journal of Sound and Vibration, 257(4), (2002), 605-613.
8. O’Connor W, Lang D, Position control of flexible robot arms using mechanical waves, Transactions of the ASME, 120 (1998) 334- 339.
9. Karagülle H, Malgaca L, Öktem H F, Analysis of active vibration control in smart structures by ANSYS., Smart Mater. Struct., 13 (2004) 661-667.
