Dynamická syntéza a redukce počtu stupňů volnosti soustav Zpracoval Doc. RNDr. Zdeněk Hlaváč, CSc
Redukcí počtu stupňů volnosti (též kondenzací) matematického modelu mechanické soustavy rozumíme jeho transformaci z prostoru dimenze n do prostoru dimenze m < n, při (přibližném) zachování vlastních frekvencí zájmového rozsahu. Velké procento metod redukce se opírá o transformaci vektoru q zobecněných souřadnic na redukovaný vektor x podle vztahu q = Tx,
(1)
kde T je na čase nezávislá matice typu (n, m). Matematický model M q¨ + B q˙ + Kq = f
(2)
potom přejde do tvaru ¨ + BT q˙ + KT x = f . MT x Násobením zleva maticí T T dostáváme redukovaný model v prostoru souřadnic x ve tvaru
kde
e, gx fx f =f ¨ +B ˙ + Kx M
(3)
e = TTf . f = T T XT , X = M , B, K , f X
(4)
Jednotlivé konkrétní metody se liší volbou transformační matice T . Úlohou dynamické syntézy rozumíme konstrukci matematického (popřípadě i dynamického) modelu soustavy pro zadané vlastní hodnoty. Jedná se tedy o inverzní úlohu k úloze modální analýzy soustavy, ve které na základě znalosti matematického modelu určujeme vlastní hodnoty. Metody dynamické syntézy mohou sloužit též k redukci počtu stupňů volnosti soustav, když pro jedenkrát modálně analyzovanou soustavu vybereme pouze vlastní čísla (a jim odpovídající vlastní vektory) pouze ze zájmového rozsahu. Koeficientové matice matematického modelu se v těchto případech konstruují ”od samého počátku” bez použití transformačních matic na matice plného modelu. Z každé skupiny metod uvedeme podrobněji dvě metody aplikované na slabě nekonzervativní soustavy. V závěru kapitoly uvedeme ještě speciální (hojně využívanou) metodu redukce počtu stupňů volnosti soustavy složené z podsoustav vázaných pružněviskózními vazbami. Tato metoda nespadá přímo do žádné ze dvou popsaných skupin metod. 1
1
Modální redukce
Nechť model (2) je slabě nekonzervativní. Proveďme modální analýzu k němu přidruženého konzervativního modelu, tedy modelu M q¨ + Kq = 0 . Vzestupně řazené kvadráty vlastních frekvencí Ω2i uspořádáme do diagonální spektrální matice Λ a jim odpovídající vlastní vektory v i po sloupcích do modální matice V . Jestliže vlastní vektory normalizujeme M -normou, platí pro modální matici diagonalizační vztah V TMV = E ,
(5)
kde E je jednotková matice takového řádu, kolik má soustava stupňů volnosti. Z definičního vztahu pro problém vlastních hodnot pak plyne (viz příslušné téma v předmětu Matematická teorie kmitání) pro modální matici též druhý diagonalizační vztah V T KV = Λ .
(6)
Vzhledem ke slabé nekonzervativnosti soustavy platí navíc V T BV = diag[2Di Ωi ] .
(7)
Provedeme výběr zájmového rozsahu vlastních frekvencí Ωp+1 , . . . , Ωp+m . K nim přiřazené vlastní vektory uspořádáme po sloupcích do redukované modální matice V r = [v p+1 , . . . , v p+m ]
(8)
typu (n, m) a provedeme tzv. neúplnou modální transformaci souřadnic q = V rx .
(9)
Je zřejmé, že matice V r má charakter transformační matice T v (1). Protože pro celou modální matici platí diagonalizační vztahy (5), (6) a (7), platí analogické vztahy i pro matici redukovanou. Je tedy
V Tr M V r = E m , V Tr BV r = diag[2Dp+i Ωp+i ], V Tr KV r = diag[Ω2p+i ], i = 1, . . . , m . Redukovaný model (3) jest potom tvaru ¨ + diag[2Dp+i Ωp+i ]x˙ + diag[Ω2p+i ]x = V Tr f . Emx
(10)
e = Přidružený konzervativní model k (10) zřejmě (přesně) vykazuje vlastní frekvence Ω i = Ωp+i , i = 1, . . . , m. Jeho (Eukleidovou normou normalizované) vlastní vektory xi (dimenze m) se neúplnou modální transformací (9) transformují na příslušné vlastní vektory q i ”velké” dimenze n. Platí tedy vztah q i = V r xi , i = 1, . . . , m.
2
2
Guyanova redukce
Zatímco předchozí metoda vyžadovala modální analýzu původní soustavy o velkém počtu stupňů volnosti, právě popisovaná metoda tento krok nutně nevyžaduje. Místo zájmového rozsahu zachovávaných vlastních frekvencí metoda vyžaduje výběr m < n tzv. hlavních souřadnic q m . Tyto souřadnice vybíráme zejména s ohledem na následující skutečnosti: • popisují pohyb dominantní hmotnosti, • jsou buzené • vyskytují se ve výpočtu deformace vazeb, jimiž přenášené silové účinky chceme finálně zjišťovat. Nyní permutujeme řádky i sloupce koeficientových matic modelu (2) tak, aby hlavní souřadnice byly na prvních m místech. Stejným způsobem se pak permutuje i vektor pravých stran. Protože do subvektoru hlavních souřadnic byly zařazeny všechny buzené souřadnice, jsou na posledních n − m místech permutovaného vektoru pravé strany nuly. Poznámka: Definujme permutační matici P ab (řádu n) jako
P ab
= En +
··· · · · −1 · · · 1 ··· ··· ··· 1 · · · −1 · · · ···
,
kde obsazená místa ve druhém sčítanci mají pořadí řádku a sloupce pouze a nebo b (a < b). Snadno se přesvědčíme, že násobení čtvercové matice řádu n zleva maticí P ab zaměňuje pořadí a−tého a b−tého řádku, zatímco násobení touto maticí zprava činí tutéž operaci se sloupci. Operací X = P k · · · P 1 XP 1 · · · P k , kde P i je výše definovaná permutační matice pro vhodně volená přirozená čísla a a b a X = M , B, K, docílíme matematického modelu tvaru M q¨ + B q˙ + Kq = f = P k · · · P 1 f ,
(11)
T
kde q T = [q Tm , q Ts ], f = [f Tm , 0T ] a koeficientové matice lze psát v blokovém tvaru X=
"
#
X mm X ms X sm X ss
, X = M , B, K .
Bloky s indexem m odpovídají hlavním souřadnicím a bloky s indexem s ostatním (tzv. vedlejším souřadnicím). Blokovým rozpisem modelu (11) dostaneme "
M mm M ms M sm M ss
#"
# "
q¨ m B B ms + mm ¨ qs B sm B ss
#"
# "
q˙ m K K ms + mm ˙q s K sm K ss
#"
# "
#
qm f = m . qs 0
(12)
Poznámka: Protože (samozřejmě) předpokládáme symetrické koeficientové matice mo3
delu, platí zřejmě X sm = X Tms , X = M , B, K . Podstatou metody je fakt, že druhou blokovou řádku modelu (12) bereme kvazistaticky, tedy zanedbáváme v ní setrvačné a tlumící členy. Tato řádka má v tom případě jednoduchý tvar K sm q m + K ss q s = 0 .
(13)
Za předpokladu regularity (čtvercové) matice K ss odtud plyne q s = −K −1 ss K sm q m , takže q=
"
qm qs
#
=
"
E mqm −K −1 ss K sm q m
#
=
"
Em −K −1 ss K sm
#
(14)
qm .
Získali jsme tak vztah mezi původními a hlavními souřadnicemi vztah tvaru (1), kde T =
"
Em −K −1 ss K sm
#
(15)
.
Podle (4) pak pro redukované koeficientové matice platí T
X r = T XT = [E m ,
K ms K −1 ss ]
"
X mm X ms X sm X ss
#"
Em −K −1 ss K sm
#
=
(16)
−1 −1 −1 = X mm − X ms K −1 ss K sm − K ms K ss X sm + K ms K ss X ss K ss K sm , X = M , B, K .
Pro X = K se výpočet zřejmě zjednoduší na tvar K r = K mm − K ms K −1 ss K sm .
(17)
Poznámka: 1. Druhou blokovou řádku (12) lze brát ve tvaru (13) v případě, že blok M mm je oproti blokům M ms a M ss dominantní a zároveň i blok B mm jest dominantní vůči blokům B ms a B ss . Toho se obecně dosahuje volbou těch zobecněných souřadnic do hlavních, kolem kterých je soustředěno tlumení a (jak už bylo výše řečeno) které popisují pohyb velkých hmotností. Tímto postupem se automaticky uvažuje (přibližně) zachovaný zájmový rozsah m nejnižších vlastních frekvencí. 2. Roznásobené výrazy (16) a (17) pro redukované matice se obvykle v praxi nepoužívají. Většinou se konkrétně napočítá matice T a z ní pak redukované koeficientové matice. 3. Jestliže soustava je bez tlumení s diagonálně dominantními maticemi M a K, postačuje vypočítat poměry mkiiii a do vektoru hlavních souřadnic zařadíme všechny souřadnice pořadí i, pro něž poměr nepřevýší vhodně zvolenou hranici. Kdyby totiž zmíněné matice byly diagonální, představovaly by tyto poměry přímo vlastní čísla úlohy, tedy kvadráty vlastních frekvencí soustavy. 4
Pokud vedlejší souřadnice bude pouze jediná, je
K ss = [knn ] ⇒ K −1 ss =
1 1 ⇒ K ms K −1 ss = knn knn
k1n .. .
.
kn−1,n
Označíme-li tento vektor jako t, obdržíme podle (16) a (17)
X r = X mm − tX Tms − X ms tT + tX ss tT , X = M .B , K r = K mm − tK Tms . Poznámka: Zřejmě platí X ss = xnn (jediný prvek), X sm = [xn1 , . . . , xn,n−1 ], X ms = = X Tsm , X = M , B, x = m, b. Proto výrazy tX Tms , X ms tT i tX ss tT = xnn ttT jsou násobením sloupcového s řádkovým vektorem. Jedná se o tzv, dyadický součin, tedy matici (řádu n − 1). Rozepíšeme-li předchozí vztahy po prvcích, dostáváme pro prvky redukovaných matic (r) (r) X r = [xij ], X = M , B, x = m, b, resp. K r = [kij ] vztahy (r)
xij = xij −
kin xjn xin kjn kin xnn kjn − + , x = m, b , 2 knn knn knn
(18)
kin kjn , i, j = 1, . . . , n − 1 . knn
(19)
(r)
kij = kij −
Procesů (18), (19), které připomínají Gaussovu eliminaci při řešení soustavy lineárních algebraických rovnic, provedeme celkem n − m. Příklad : Mějme přímý zleva vetknutý řetězec se třemi stupni volnosti s hmotovými parametry mi a tuhostními parametry ki , i = 1, 2, 3 podle obr.1. K němu je prostřednictvím tuhosti k0 připojena řádově menší hmotnost m0 . Redukujeme vzniklou soustavu se čtyřmi stupni volnosti o jeden stupeň Guyanovou redukcí zavedením vedlejší souřadnice popisující pohyb ”malé” hmotnosti.
k1
m1
k2
m2
k3
m3
k0
m0
Obrázek 1: Řešení: Matematický model soustavy se čtyřmi stupni volnosti jest zřejmě popsán maticemi hmotnosti M a tuhosti K tvaru
M = diag[m1 , m2 , m3 , m0 ] , K =
k1 + k2 −k2 0 0 −k2 k2 + k3 −k3 0 0 −k3 k3 + k0 −k0 0 0 −k0 k0
.
−1 1 Dále platí K ss = [k0 ], takže K −1 ss = [ k0 ]. Dále jest K sm = [0, 0, −k0 ], pročež −K ss K sm =
= [0, 0, 1]. Transformační matice T má podle (15) tvar 5
T =
1 0 0 0
0 1 0 0
0 0 1 1
Pro matice redukovaného modelu potom platí
1 0 0 0 Mr = T TMT = 0 1 0 0 0 0 1 1
m1 0 0 0 0 = 0 m2 0 0 0 m3 m0
1 0 0 0 0 0 0 0 1 1
K r = T KT = 0 1 T
1 0 0 0
.
m1 0 0 0 0 m2 0 0 0 0 m3 0 0 0 0 m0
0 1 0 0
0 0 1 1
m1 = 0
0
1 0 0 0
0 1 0 0
0 0 1 1
0 1 0 0
0 0 1 1
=
0 0 m2 0 , 0 m3 + m0
k1 + k2 −k2 0 0 −k2 k2 + k3 −k3 0 0 −k3 k3 + k0 −k0 0 0 −k0 k0
1 k1 + k2 −k2 0 0 0 k2 + k3 −k3 0 = −k2 0 0 −k3 k3 0 0
k1 + k2 = −k2
0
1 0 0 0
0 1 0 0
0 0 1 1
=
−k2 0 k2 + k3 −k3 . −k3 k3
Vzniklý redukovaný model se třemi stupni volnosti (obr.2) zřejmě odpovídá dynamickému modelu zleva vetknutého přímého řetězce s hmotovými parametry m1 , m2 a m3 + m0 a tuhostními parametry k1 , k2 , k3 . Na tuhosti k0 redukovaný model vůbec nezávisí. ”Malá” hmotnost m0 je ”natvrdo” připojena ke hmotnosti m3 . m3+m k1
m1
k2
m2
0
k3
Obrázek 2: Představu o chybě náhrady redukovaným modelem si lze udělat z číselného zadání m1 = 1, m2 = 2, m3 = 3, m0 = 0.1 [kg] a k1 = 104 , k2 = 2 · 104 , k3 = 5 · 104 , k0 = 105 [N/m]. Modální analýzou získáme vlastní frekvence původní soustavy jako Ω1 = 33.903, Ω2 = 161.20, Ω3 = 232.32, Ω4 = 1016.8 [rad/s]. Vlastní frekvence redukované soustavy pak jsou Ω1 = 33.903, Ω2 = 161.22, Ω3 = 232.35 [rad/s], tudíž aproximují tři nejnižší vlastní frekvence původní soustavy s uspokojivou přesností. Zadáme-li m0 = 0.01 [kg] (a ostatní parametry stejné), vyjde Ω1 = 34.198, Ω2 = 161.78, Ω3 = 232.96, Ω4 = 3167.6 [rad/s] a Ω1 = 34.198, Ω2 = 161.78, Ω3 = 232.96 [rad/s]. Aproximace tří nejnižších vlastních frekvencí původní soustavy je nyní ještě lepší. Poznámka: Analogický výsledek bychom (technicky podstatně náročnější cestou) dostali i pro řetězec s obecně n stupni volnosti. Provedeme-li na redukované soustavě další 6
redukci, dostáváme řetězce o stále menším počtu stupňů volnosti se stejnými parametry, pouze k poslední hmotnosti se přičítají hmotnosti charakterizující změny všech vedlejších souřadnic. Jestliže tedy máme přímý řetězec s n stupni volnosti a posledních n − m hmotností je řádově nižších než prvních m hmotností (obr.3), dospějeme Guyanovou redukcí postupně až k modelu na obr.4, s parametry (r)
(r)
mj = mj , j = 1, . . . , m − 1 , kj = kj , j = 1, . . . , m , m(r) m =
n X
mi .
i=m
m m+1 k1
km
m1
k m+1
mm
kn
mn
Obrázek 3:
n
mi
m m−1 k1
km
k m−1
m1
i=m
Obrázek 4:
3
Metoda inverzních vzorců (matematické kondenzace)
Mějme slabě nekonzervativní model (2). Provedením modální analýzy přidruženého konzervativního modelu získáme spektrální matici Λ = diag[Ω2i ] složenou z kvadrátů vlastních frekvencí (přidruženého konzervativního modelu) a modální matici V složenou z M −normou normalizovaných vlastních vektorů k těmto frekvencím přiřazených. Pak zřejmě platí V T M V = E ; V T KV = Λ = diag[Ω2i ] ; V T BV = D = diag[2Di Ωi ] ,
(20)
kde Di je poměrný útlum příslušející k i−tému tvaru kmitání. Protože vlastní vektory tvoří bázi Eukleidova prostoru Rn , existuje inverzní modální matice V −1 a tudíž i inverzní k transponované modální matici (V T )−1 stručně označovaná V −T . Násobením rovnic (20) zleva maticí V −T a zprava maticí V −1 dostaneme M = V −T V −1 = (V V T )−1 ; B = V −T DV −1 = (V D −1 V T )−1 ; K=V diag[ 2D1i Ωi ]
−T
ΛV
−1
T −1
= (V Λ V ) −1
diag[ Ω12 ]. i
(21)
,
aΛ = Výrazům (21) říkáme inverzní vzorce. Řeší kde D = přesně úlohu dynamické syntézy, tedy na základě zadaných vlastních hodnot přidružené konzervativní soustavy a poměrných útlumů určit její matematický model. Metoda −1
−1
7
má ovšem jednu podstatnou nevýhodu. Koeficientové matice získáme plné, takže jejich jednotlivým prvkům nelze přiřadit hmotový parametr ani tuhostní a tlumící parametr vazby. Ke vzniklému matematickému modelu tedy nelze uspokojivě přiřadit model dynamický. Modifikované inverzní vzorce lze využít též pro redukci počtu stupňů volnosti už zadaného dynamického (matematického) modelu. Vybereme zájmovou oblast frekvencí, do níž spadají vlastní frekvence Ωp+1 , . . . , Ωp+m modálně analyzovaného přidruženého konzervativního plného modelu. Do matice V umístíme pouze vlastní vektory přiřazené vybraným frekvencím, tedy pořadových čísel p + 1, . . . , p + m. Matice V je proto typu (n, m). Její prvky označme v ij . Pro každou řádku této matice utvořme součet absolutních hodnot jejích prvků µi =
m X
|v ij |, i = 1, . . . , n. Za vybrané řádky (místa kondenzace)
j=1
volíme m řádků, kde parametry µi nabývají m největších hodnot. Nechť jsou to řádky pořadových čísel i1 , . . . , im . Pak definujeme spektrální a modální matici redukovaného modelu jako Λr = diag[Ωp+k ] ; V r = [v ik ,j ], j, k = 1, . . . , m .
(22)
Pro zadané poměrné útlumy Dp+1 , . . . , Dp+m pak určíme matici D r = diag[2Dp+i Ωp+i ] a podle modifikovaných inverzních vzorců T −1 T −1 M r = (V r V Tr )−1 ; B r = (V r D −1 ; K r = (V r Λ−1 r Vr) r Vr)
(23)
vyčíslíme koeficientové matice redukovaného modelu s m stupni volnosti. Příklad : Mějme soustavu tvaru zleva vetknutého přímého řetězce s deseti stupni volnosti o parametrech mi = i [kg] a ki = i · 104 [N/m] (i = 1, . . . , 10). Redukujme tuto soustavu pomocí inverzních vzorců na soustavu se třemi stupni volnosti, při zachování nejnižších tří vlastních frekvencí původního plného (netlumeného) modelu. Řešení: Matematický model soustavy má diagonální matici hmotnosti M = diag[mi ] a tridiagonální (symetrickou) matici tuhosti K, jež má na diagonále vektor dT = = [k1 +k2 , k2 +k3 , . . . , k9 +k10 , k10 ] a na přidiagonálách vektor pT = [−k2 , −k3 , . . . , −k10 ]. Modální analýzou tohoto modelu získáme vlastní frekvence soustavy jako Ω1 = 8.807, Ω2 = 44.49, Ω3 = 76.40, Ω4 = 106.2, Ω5 = 133.3, Ω6 = 157.2, Ω7 = 177.2, Ω8 = 192.9, Ω9 = 203.5, Ω10 = 216.0 [rad/s] . Neúplná modální matice vlastních vektorů přiřazených třem nejnižším vlastním frekvencím (normalizovaná M − normou) je
Vr
=
−5.6982e − 02 1.6377e − 01 −8.5252e − 02 2.2945e − 01 −1.0366e − 01 2.4295e − 01 −1.1686e − 01 2.1701e − 01 −1.2669e − 01 1.6189e − 01 −1.3407e − 01 8.9249e − 02 −1.3950e − 01 1.1844e − 02 −1.4331e − 01 −5.7938e − 02 −1.4571e − 01 −1.0977e − 01 −1.4685e − 01 −1.3686e − 01 8
−2.4343e − 01 −2.9411e − 01 −2.1346e − 01 −5.9537e − 02 9.1398e − 02 1.7273e − 01 1.5603e − 01 6.1738e − 02 −5.4104e − 02 −1.2994e − 01
.
Vybereme nyní tři souřadnice, jež při kmitání prvními třemi tvary kmitání ”nejvíce kmitají”. Utvoříme tedy vektor součtů absolutních hodnot prvků matice V r po řádcích. Dostaneme vektor sT = [0.464, 0.609, 0.560, 0.393, 0.380, 0.396, 0, 307, 0.263, 0.310, 0.414] . Srovnáním prvků tohoto vektoru podle velikosti vidíme, že ”nejvíc kmitá” při kmitání prvními třemi tvary kmitání druhá, třetí a první souřadnice. Tyto souřadnice tedy bereme jako ”místa kondenzace”, takže redukovaná modální matice V r a redukovaná spektrální matice Λr mají tvar
−5.6982e − 02 1.6377e − 01 −2.4343e − 01 Vr = −8.5252e − 02 2.2945e − 01 −2.9411e − 01 ; −1.0366e − 01 2.4295e − 01 −2.1346e − 01 Λr = diag[Ω21 , Ω22 , Ω23 ] = diag[77.561, 1979.6, 5836.3] .
Dosazením do modifikovaných inverzních vzorců (23) získáme redukovanou matici hmotnosti M r a tuhosti K r tvaru
2.6021e + 06 −2.7903e + 06 8.7801e + 05 −2.7903e + 06 2.9923e + 06 −9.4171e + 05 Mr = ; 8.7801e + 05 −9.4171e + 05 2.9648e + 05 1.5371e + 09 −1.6272e + 09 4.9440e + 08 Kr = −1.6272e + 09 1.7231e + 09 −5.2375e + 08 . 4.9440e + 08 −5.2375e + 08 1.5933e + 08
Je patrno, že obě matice redukovaného modelu jsou plné, takže nelze k nim přiřadit vhodný dynamický model. Poznámka: Provedením modální analýzy získaného redukovaného modelu obdržíme přesně první tři vlastní frekvence původního modelu s deseti stupni volnosti i vlastní vektory skládající se z prvních tří souřadnic prvních tří vlastních vektorů původního modelu.
4
Metoda redukce na model předepsané struktury (fyzikální kondenzace)
Předpokládejme opět slabě tlumenou soustavu s velkým počtem n stupňů volnosti. Provedeme modální analýzu přidružené konzervativní soustavy. Její vlastní frekvence nechť jsou Ω∗i a j−tá souřadnice i−tého vlastního vektoru (normalizovaného M −normou) nechť je vji , i, j = 1, . . . , n. Vybereme zájmový interval (přibližně) zachovávaných vlastních frekvencí Ωr+1 , . . . , Ωr+p , p << n. Vyberme dále tzv. místa kondenzace, tedy zobecněné souřadnice qij , j = 1, . . . , m, m << n, p ≤ m, jež budou odpovídat souřadnicím popisujícím konfiguraci redukovaného modelu dané topologické struktury. Jejich výběr bude určen vektorem výběru ij , j = 1, . . . , m. Poznámka: Obvykle volíme za zájmový interval zachovávaných vlastních frekvencí (pro neizolované soustavy) p nejnižších vlastních frekvencí, tedy volíme r = 0. Místa kondenzace volíme s ohledem na zachování základních rysů topologické struktury původního 9
(velkého) modelu. Volba míst kondenzace jako souřadnic popisujících pohyb velkých hmotností (jako u Guyanovy redukce) se zde nejeví příliš vhodná. Více o volbě míst kondenzace viz příklad níže. Redukovaný model bude vykazovat m vlastních frekvencí Ω∗i = Ωr+i , i = 1, . . . , m a ∗ jeho vlastní vektory v ∗i budou mít souřadnice vj,i = vij ,r+i , i, j = 1, . . . , m. Hvězdičky v označení vlastních hodnot redukovaného modelu budeme v dalším vynechávat, protože vlastní hodnoty velkého modelu (bez hvězdiček) už nebudeme potřebovat. Seřaďme hmotové parametry redukovaného modelu do vektoru hmotových parametrů m dimenze sm a tuhostní parametry do vektoru tuhostních parametrů k dimenze sk . Pomocí těchto parametrů jsou sestaveny matice hmotnosti redukovaného modelu M R a jeho matice tuhosti K R . Budeme požadovat, aby vlastní vektory redukovaného modelu byly normovány M R −normou, tedy aby platilo v Tk M R v l = δkl ; v Tk K R v l = Ω2k δkl
(24)
kde δkl je Kronekerův symbol. Protože prvky matice M R závisejí lineárně na prvcích vektoru hmotnosti m, a prvky matice K R zase závisejí lineárně na prvcích vektoru tuhosti k, existují pro každé l ∈ {1, . . . , m} matice X l a matice Y l (složené ze souřadnic vlastních vektorů), že platí M Rvl = X l m ; K Rvl = Y l k .
(25)
Matice X l je zřejmě typu (m, sm ) a matice Y pak typu (m, sk ). Diagonalizační vztahy (24) lze pomocí (25) přepsat do tvaru v Tk X l m = δkl ; v Tk Y l k = Ω2k δkl ,
(26)
pro k, l = 1, . . . , m. Vzhledem k symetrii matic M R a K R stačí (26) zapsat pouze pro taková pořadí řádků a sloupců, jež vyjadřují horní trojúhelník výsledných matic (včetně hlavní diagonály) pro p nejnižších vlastních frekvencí redukované soustavy. Výrazy (26) stačí tedy např. psát pro k = 1, . . . , p a l = k, . . . , p. Zapíšeme-li soustavy rovnic (26) pro různé dvojice možných indexů k a l (např. tak, že sloupcový index l se mění rychleji), dostáváme po vytknutí z prvních soustav vektoru hmotnosti m a ze druhých soustav vektoru tuhosti k soustavy rovnic tvaru
v T1 X 1 v1X 2 .. . v T1 X p v2X 2 v2X 3 .. . v2X p .. . vpX p
1 0 .. .
v T1 Y 1 v1Y 2 .. .
0 vT Y p 1 1 v2Y 2 m = 0 ; v2Y 3 . .. .. . 0 v Y 2 p .. .. . .
1
vpY p
Ω21 0 .. .
0 2 Ω k = 2 0 . .. 0 .. .
Ω2p
.
(27)
Každý blok v maticích na levých stranách je typu (1, sm ) resp. (1, sk ) a bloků je celkem p (p+1), tedy tolik, kolik jest různých dvojic v úvahu připadajících indexů k a l. Rovnice 2 přepíšeme symbolicky do tvaru Φm = δ 1 ; Ψk = δ 2 . 10
(28)
Matice Φ je obecně obdélníková typu ( p2 (p + 1), sm ) a matice Ψ je rovněž obecně obdélníková typu ( p2 (p + 1), sk ). Soustavy (28) jsou tedy soustavami lineárních algebraických rovnic pro neznámé vektory hmotnosti m a tuhosti k redukovaného modelu dané topologické struktury popsané maticemi M R a K R . Kvalita soustav rovnic závisí na volbě počtu p zachovávaných vlastních frekvencí. Pokud volíme p = m jsou soustavy (28) ve většině případů přeurčené, neboť je zpravidla m2 (m + 1) > max(sm , sk ). Jejich přesné řešení pak ve velké většině případů neexistuje (viz příslušné téma). Formulovat lze pouze ”nejlepší” pseudořešení ve tvaru m = ΦL δ 1 ; k = ΨL δ 2 ,
(29)
ΦL = (ΦT Φ)−1 ΦT ; ΨL = (ΨT Ψ)−1 ΨT
(30)
kde
jsou tzv. levé pseudoinverzní matice (viz příslušné téma). Tyto matice existují za předpokladu regularity matic Φ a Ψ, kdy jejich sloupce jsou lineárně nezávislé. Mají tu vlastnost, že utvoříme-li podle (29) vektory m resp. k, pak ||Φm − δ 1 || = min ||Φx − δ 1 || ; ||Ψk − δ 2 || = min ||Ψx − δ 2 || . x∈Rsm x∈Rsk
(31)
Eukleidova norma reziduálního vektoru, tedy odchylky dosazeného pseudořešení (29) do levých stran rovnic (28) od odpovídajících pravých stran těchto rovnic je nejmenší ze všech možných dosazení vektorů odpovídající dimenze. Volíme-li p < m (asi ve většině praktických případů), mohou být soustavy (28) i nedourčené, kdy počet rovnic p2 (p + 1) je menší než počet neznámých sm resp. sk . V tom případě zpravidla řešení soustav (28) je nekonečná množina. Vyjádřeme jedno z nich ve tvaru m = ΦP δ 1 ; k = ΨP δ 2 ,
(32)
ΦP = ΦT (ΦΦT )−1 ; ΨP = ΨT (ΨΨT )−1
(33)
kde
jsou tzv. pravé pseudoinverzní matice (viz příslušné téma). Tyto matice existují za předpokladu regularity matic Φ a Ψ, kdy jejich řádky jsou lineárně nezávislé. Řešení (32) má potom ze všech v úvahu připadajících řešení soustav (28) nejmenší Eukleidovu normu. Ve výjimečných případech se může stát, že p2 (p + 1) = sm nebo p2 (p + 1) = sk . V takovém případě pro regulární matici Φ resp. Ψ existuje jediné řešení soustav (28) ve tvaru m = Φ−1 δ 1 ; m = Ψ−1 δ 2 ,
(34)
protože příslušné inverzní matice existují. Speciálně pro případ, že redukovaná soustava má topologickou strukturu torzní (resp. podélně kmitající) soustavy popsané incidenční maticí I = [ijk ], j = 1, 2; k = 1, . . . , sk o hmotových parametrech uspořádaných do vektoru m = [m1 , . . . , msm ]T a tuhostních parametrech uspořádaných do vektoru k = [k1 , . . . , ksk ]T , lze matice Φ a Ψ určit přímo prostřednictvím souřadnic vlastních vektorů redukovaného modelu (výše označované hvězdičkou). Pak totiž matice matematického modelu redukované soustavy má tvar 11
M R = diag[m1 , . . . , msm ] ; K R
= kj j=1 sk X
... . . . 1 . . . −1 . . . ... . . . −1 . . . 1 . . . ...
.
(35)
Pozice zaplněné nenulovými prvky ve sčítancích matice tuhosti tvoří hnízdo o pořadových číslech řádků a sloupců i1,j a i2,j . Pokud i1j = 0, jest zaplněn jedničkou pouze diagonální prvek na pozici (i2,j , i2j ). Z prvního výrazu (35) plyne, že pro libovolné l je (hvězdičky u vlastních vektorů a jejich souřadnic pro redukovaný model vynecháváme a uvažujeme modální matici už pouze typu (sm , p), kdy její řádky a sloupce jsou vybrané podle předpisu výše)
m1 v1l v1l . .. = diag[v1l , . . . , vs l ]m . M R v l = diag[m1 , . . . , msm ] . m .. = m sm v sm l v sm l Odtud plyne, že X l = diag[v1l , . . . , vsm l ] . Navíc pro libovolné k platí v Tk M R v l = [v1k , . . . , vsm k ]diag[v1l , . . . , vsm l ]m = [v1k v1l , . . . , vsm k vsm l ]m . Podle definice matice Φ tato má zřejmě tvar
Φ=
2 v11 v11 v12 ... v11 v1p 2 v12 ... v12 v1p ... 2 v1p
... vs2m 1 . . . v sm 1 v sm 2 . . . v sm 1 v sm p ... vs2m 2
. . . v sm 2 v sm p vs2m p
...
Z druhého výrazu (35) plyne, že pro libovolné l je
K Rvl = kj j=1 sk X
... . . . 1 . . . −1 . . . ... . . . −1 . . . 1 . . . ...
.
(36)
v1l vi1j l − vi2j l sk X . . = k j . j=1 vi2j l − vi1j l v sm l
,
přičemž zaplněná místa mají pořadová čísla (řádků) i1j a i2j , j = 1, . . . , sk . Po sečtení přes všechny tuhostní vazby redukovaného modelu dostáváme zřejmě sloupec o sm prvcích výše uvedeného charakteru. Každý prvek odpovídá hmotě a obsahuje součty sčítanců výše uvedeného charakteru, z nichž každý odpovídá tuhostní vazbě, která inciduje s danou hmotou. Odtud vyplývá, že sloupec K R v l lze zapsat ve tvaru 12
vi11 l − vi21 l 0 K Rvl = v i21 l − vi11 l
...
0
. . . vi1sk l − vi2sk l ...
0
0
. . . vi2sk l − vi1sk l
k.
Ve zde rozepsané matici má j−tý sloupec (j = 1, . . . , sk ) zaplněna místa pořadových čísel i1j a i2j . Pro libovolné k potom platí
vi11 l − vi21 l 0 T v k K R v l = [v1k , . . . , vsm k ] v i21 l − vi11 l
0
...
0
. . . vi1sk l − vi2sk l ...
0
. . . vi2sk l − vi1sk l
k =
= [(vi11 k − vi21 k )(vi11 l − vi21 l ), . . . , (vi2sk k − vi1sk k )(vi2sk l − vi1sk l )]k . Podle definice matice Ψ tato má zřejmě tvar
Ψ=
(vi11 1 − vi21 1 )2 (vi11 1 − vi21 1 )(vi11 2 − vi21 2 ) ... (vi11 1 − vi21 1 )(vi11 p − vi21 p ) ... (vi11 2 − vi21 2 )2 ... (vi11 2 − vi21 2 )(vi11 p − vi21 p ) ... (vi11 p − vi21 p )2
... (vi2sk 1 − vi1sk 1 )2 . . . (vi2sk 1 − vi1sk 1 )(vi2sk 2 − vi1sk 2 ) . . . (vi2sk 1
... . . . (vi2sk 2 ...
− vi1sk 1 )(vi2sk p − vi1sk p ) . 2 (vi2sk 2 − vi1sk 2 ) − vi1sk 2 )(vi2sk p − vi1sk p )
(37)
(vi2sk p − vi1sk p )2
Ukážeme nyní dva konkrétní příklady redukce počtu stupňů volnosti naznačenou metodikou. První příklad ukáže pouze soustavy rovnic pro výpočet hmotových a tuhostních parametrů redukované soustavy bez konkrétních číselných výsledků. Druhý příklad ukáže na jednoduché soustavě prostřednictvím číselně zadaných parametrů úskalí volby míst kondenzace. Příklad 1: Mějme rozvětvený (např. torzní) řetězec s devíti stupni volnosti, topologické struktury podle obr.5. V místech kondenzace I2 , I4 a I6 proveďme jeho fyzikální kondenzaci na model analogické topologické struktury, znázorněný na obr.6, tak, aby byly (přibližně) zachovány dvě nejnižší vlastní frekvence původního modelu. Je tedy v tomto příkladu n = 9, m = 3, r = 0, p = 2 a i = [ij ] = [2, 4, 6]T . Řešení: Provedeme modální analýzu původní soustavy, čímž dostaneme vlastní frekvence Ωi a (příslušným způsobem normované) vlastní vektory v i = [vji ], i, j = 1, . . . , 9. 13
k1
I1
k2
k3
I3
k6
I6
k4
k5
I4
I5
I2 k7
k8
I7
I8
k9
I9
Obrázek 5:
+
k1
+
k2
+
I2
k +3
I3
I1 +
Obrázek 6: Volíme pro kondenzovaný model zachované vlastní frekvence Ω∗i = Ωi , i = 1, 2 a pří∗ ∗ ∗ slušné vlastní vektory v ∗i = [vki ], i = 1, 2, k = 1, 2, 3, kdy v1i = v2i , v2i = v4i a ∗ v3i = v6i , i = 1, 2. Vzhledem k zadané topologické struktuře kondenzovaného modelu o hmotových parametrech Ii∗ , i = 1, 2, 3 a tuhostních parametrech ki∗ , i = 1, 2, 3 má matice hmotnosti M R a matice tuhosti K R tohoto modelu tvar
I1∗ 0 0 k1∗ + k2∗ + k3∗ −k2∗ −k3∗ −k2∗ k2∗ 0 MR = 0 I2∗ 0 ; K R = . ∗ ∗ ∗ 0 0 I3 −k3 0 k3
V dalším textu už veličiny původního modelu nebudeme potřebovat, takže hvězdičky v popisu veličin kondenzovaného modelu budeme vynechávat. Pro libovolné l = 1, 2, 3 tak máme M R v l = [I1 v1l , I2 v2l , I3 v3l ]T , K R v l = [(k1 + k2 + k3 )v1l − k2 v2l − k3 v3l , −k2 v1l + k2 v2l , −k3 v1l + k3 v3l ]T = = [k1 v1l + k2 (v1l − v2l ) + k3 (v1l − v3l ), k2 (v2l − v1l ), k3 (v3l − v1l )]T . S ohledem na definici matic X l a Y l v (25) odtud máme pro ně výrazy
v1l 0 0 M R v l = X l m = 0 v2l 0 m, 0 0 v3l
v1l v1,l − v2l v1l − v3l 0 K Rvl = Y l k = 0 v2l − v1l k. 0 0 v3l − v1l
Násobením těchto vztahů zleva vektorem v k odtud vzhledem k předpokládané platnosti diagonalizačních vztahů dostáváme postupně 14
v1l 0 0 T δkl = v k M R v l = v k X l m = [v1k , v2k , v3k ] 0 v2l 0 m = 0 0 v3l = [v1k v1l , v2k v2l , v3k v3l ]m
a analogicky
v1l v1,l − v2l v1l − v3l 0 δkl Ω2k = v Tk K R v l = v k Y l k = [v1k , v2k , v3k ] k = 0 v2l − v1l 0 0 v3l − v1l = [v1k v1l , v1k (v1l − v2l ) + v2k (v2l − v1l ), v1k (v1l − v3l ) + v3k (v3l − v1l )]k = = [v1k v1l , (v1k − v2k )(v1l − v2l ), (v1k − v3k )(v1l − v3l )]k . Sepíšeme-li první a poslední vztah předchozích dvou soustav rovností pod sebe vždy pro k = 1, 2 a l = k, . . . , 2, přičemž index l se mění rychleji, dostaneme, s ohledem na definici matic a vektorů ve vztazích (28) formálně stejné rovnice pro výpočet parametrů kondenzovaného modelu, ve kterých ovšem jest
2 2 2 v11 v21 v31 Φ = v11 v12 v21, v22 v31 v32 ; 2 2 2 v12 v22 v32
2 v11 (v11 − v21 )2 (v11 − v31 )2 Ψ = v11 v12 (v11 − v21 )(v12 − v22 ) (v11 − v31 )(v12 − v32 ) ; 2 v12 (v12 − v22 )2 (v12 − v32 )2
δ T1 = [1, 0, 1] ; δ T2 = [Ω21 , 0, , Ω22 ] .
Soustavy (28) se zde určenými maticemi soustav a vektory pravých stran jsou tedy soustavami tří rovnic pro tři neznámé momenty setrvačnosti a tři neznámé torzní tuhosti kondenzovaného modelu z obr.6. Příklad 2: Mějme přímý řetězec se sedmi stupni volnosti topologické struktury podle obr.7 o hmotových parametrech m1 = m2 = 10[kg], m3 = m4 = m5 = m6 = m7 = 1[kg] a tuhostních parametrech ki = i · 104 [N/m] (i = 1, . . . , 7). Proveďme redukci soustavy na řetězec analogické topologické struktury (tedy ”zleva” vetknutý) se dvěma stupni volnosti tak, aby se (přibližně) zachovaly dvě nejnižší vlastní frekvence původní soustavy. Proveďme numerický experiment s místy kondenzace. k1
m1
k2
m2
k3
m3
k4
m4
k5
k6
m5
m6
k7
m7
Obrázek 7: Řešení: Proveďme modální analýzu původní soustavy, jež má incidenční matici (viz příslušné téma) tvaru I = [ikl ] =
"
0 1 2 3 4 5 6 1 2 3 4 5 6 7 15
#
.
Provedením modální analýzy dostáváme vlastní frekvence modelu Ω1 = 18.02, Ω2 = 55.31, Ω3 = 77.78, Ω4 = 184.4, Ω5 = 288.6, Ω6 = 374.67, Ω7 = 459.23[rad/s]. Vlastní vektory kondenzovaného modelu potom jsou v 1 = [0.16036; 0.2145; 0.227370; 0.23517; 0.23989; 0.24252; 0.24365]T , v 2 = [0.22999; −0.0067723; −0.15771; −0.25885; −0.32393; −0.36164; −0.37817]T , v 3 = [0.1462; −0.22296; −0.019426; 0.13616; 0.24416; 0.30954; 0.33882]T , v 4 = [0.0038175; −0.05918; 0.56961; 0.55697; 0.16807; −0.25127; −0.48864]T , v 5 = [0.00063856; −0.025634; 0.66853; −0.20286; −0.56206; −0.081172; 0.4276]T , v 6 = [−0.00012354; 0.0084863; −0.38287; 0.66729; −0.36606; −0.37073; 0.36874]T , v 7 = [0.00000998; −0.0010374; 0.071195; −0.24999; 0.5475; −0.71233; 0.35391]T . Protože redukovaná soustava má dva stupně volnosti a její hmotové i tuhostní parametry jsou rovněž dva, jsou obě soustavy typu (28) soustavy tří rovnic pro dvě neznámé. Volíme-li za místa kondenzace první dvě hmoty (k čemuž by nás mohlo vést soustředění hmotnosti do těchto souřadnic), mají matice výše zmíněných soustav tvar
2 2 2 v11 v21 v11 (v11 − v21 )2 Φ = v11 v12 v21 v22 ; Ψ = v11 v12 (v11 − v21 )(v12 − v22 ) 2 2 2 v12 v22 v12 (v12 − v22 )2
a vektory pravých stran jsou
δ1 = [1, 0, 1]T ; δ2 = [Ω21 , 0, Ω22 ]T . ”Nejlepším” pseudořešením soustav rovnic výše dostaneme mR = [12.814, 14.889]T [kg], kR = [1.2932, 4.1994]T · 104 [N/m] . Modální analýzou takto redukovaného modelu modelu získáváme vlastní frekvence Ω∗1 = = 20.64, Ω∗2 = 81.74, [rad/s]. První vlastní frekvence tedy aproximuje nejnižší vlastní frekvenci plného modelu s relativní přesnosti 14.5%, zatímco druhá vlastní frekvence s vlastními frekvencemi plného modelu vůbec nesouvisí. Zmíněná redukce je tedy zcela nevhodná. Volíme-li za místa kondenzace první a třetí hmotu, mají matice soustav (28) tvar
2 2 2 v11 v31 v11 (v11 − v31 )2 Φ = v11 v12 v31 v32 ; Ψ = v11 v12 (v11 − v31 )(v12 − v32 ) 2 2 2 v12 v32 v12 (v12 − v32 )2
a vektory pravých stran jsou stejné jako v předchozím případě. ”Nejlepším” pseudořešením těchto soustav jsou vektory mR = [12.73, 13.051]T [kg], kR = [1.1085, 1.6423]T · 104 [N/m] . Modální analýzou takto redukovaného modelu modelu získáváme vlastní frekvence Ω∗1 = = 18.92, Ω∗2 = 55.33, [rad/s]. První vlastní frekvence tedy aproximuje nejnižší vlastní frekvenci plného modelu s relativní přesnosti 5% a druhá vlastní frekvence aproximuje 16
druhou nejnižší vlastní frekvenci plného modelu dokonce s relativní přesností 0.036%. Je zajímavé si na tomto místě všimnout i modální matice (maticí hmotnosti normovaných) vlastních vektorů redukované soustavy. Vychází V =
"
−0.16176 −0.22889 −0.22606 0.15976
#
.
Prvky této matice velice dobře aproximují (až na znaménko) prvky submatice modální matice plného modelu složené z prvních dvou sloupců (zachované frekvence) a prvního a třetího řádku (místa kondenzace). Poznamenejme ještě, že znaménková změna není na závadu, protože vlastní vektory soustavy jsou (i při jejich normalizaci) určeny až na znaménko, jak vyplývá z diagonalizačních vztahů v Ti M v j = δij a v Ti Kv j = Ω2i δij . Redukce s místy kondenzace v prvních dvou hmotách nebyla úspěšná z důvodů (v absolutní hodnotě) malého prvku na druhém místě druhého vlastního vektoru plného modelu. Úspěšná redukce počtu stupňů volnosti touto metodikou vyžaduje jisté zkušenosti řešitele. Obecně je možno říci, že místa kondenzace je nutno volit tak, aby byly kvalitativně zachovány tvary kmitání plného modelu při kmitání zachovávanými vlastními frekvencemi. Není-li tato podmínka splněna, stane se (zejména u soustav složitější topologické struktury), že ”řešení” soustav (28) budou mít i záporné prvky. Takovým případům jednak neodpovídá reálná redukovaná soustava a rovněž vlastní čísla takových soustav jsou obecně komplexní.
5
Modální redukce systémů složených ze subsystémů
Mechanický systém nechť je složený z N subsystémů spojených mezi sebou diskrétními, pružně-viskózními vazbami (obr.8). Každý i−tý subsystém izolovaný od ostaních nechť je (při lokálním číslování zobecněných souřadnic) charakterizován maticemi hmotnosti M i , tlumení B i , tuhosti K i a vektorem vnějších budících účinků f E i (t) (i = 1, . . . , N ). Jednotlivé subsystémy nechť jsou mezi sebou spojeny M pružně-viskózními vazbami charakterizovanými tuhostí (lineární pružiny) kjv a tlumícím koeficientem (viskózního tlumiče) bvj (j = 1, . . . , M ). Poznamenejme na tomto místě, že pokud v systému existují vazby k rámu, respektují se v maticích tuhosti a tlumení jednotlivých subsystémů.
v
bj E
E
M i , Bi , Ki , f i
M i+1 , B i+1 , K i+1 , f i+1 v
kj i − ty subsystem
j − ta vazba
(i + 1) − ty subsystem
Obrázek 8: Uvolněním každého subsystému od ostatních připojíme silové účinky f Ii ostatních subsystémů prostřednictvím přerušených spojovacích vazeb. Matematický model i−tého subsystému pak bude I M i q¨ i + B i q˙ i + K i q i = f E i (t) + f i (t) ,
17
(38)
kde q i (i = 1, . . . , N ) je vektor zobecněných souřadnic popisující konfiguraci i−tého subsystému. K modelu (38) přiřaďme příslušný konzervativní model (bez tlumení a budících účinků) a proveďme jeho modální analýzu. Získáme tak spektrální matici Λi a M i -normou normovanou modální matici V i , pro které platí V Ti M i V i = E ni ; V Ti K i V i = Λi , i = 1, . . . , N ,
(39)
kde ni je počet stupňů volnosti i−tého subsystému. Tvary kmitání každého subsystému nyní rozdělme na (obvykle přiřazené k mi nejnižším vlastním frekvencím) tzv. hlavní (m) (m) tvary kterým odpovídá spektrální submatice Λi řádu mi a modální submatice V i (o) typu (ni , mi ) a ostatní tvary kmitání, kterým odpovídají příslušné submatice Λi a (o) V i doplňkových typů. Spektrální a modální matici i−tého subsystému pak lze psát v rozděleném tvaru (m)
(o)
(m)
(o)
Λi = diag[Λi , Λi ] ; V i = [V i , V i ] . Poznámka: Jestliže bychom vybrali obecnější množinu mi hlavních tvarů, provedeme permutaci sloupců (resp. řádků i sloupců) spektrální a modální matice tak, aby hlavní tvary odpovídaly prvním mi pořadovým číslům. V každém subsystému nyní provedeme neúplnou modální transformaci (m)
q i = V i xi ; i = 1, . . . , N ,
(40)
kde xi je vektor hlavních modálních souřadnic i−tého subsystému. Dosazením (40) do modelu (38) dostaneme (m) E I ˙ i + K i V (m) ¨ i + B i V (m) M iV i x i xi = f i (t) + f i (t) . i x
(41)
Protože platí (39), platí zřejmě též (m)T
Vi
(m)
M iV i
(m)T
= E mi ; V i (m)T
Násobením rovnosti (41) zleva maticí V I ních výrazů obdržíme (m)T
¨i + V i E mi x
(m)
K iV i
(m)
= Λi
.
a uvážením výše uvedených diagonalizač-
(m) (m)T (m) I B i V i x˙ i + Λi xi = V i (f E i (t) + f i (t)) .
(42)
Zaveďme nyní globální vektory celého systému jako q T = [q T1 , . . . , q TN ], xT = [xT1 , . . . , xTN ], E T T I I T I T T T f E = [(f E 1 ) , . . . , (f N ) ] , f = [(f 1 ) , . . . , (f N ) ]
a matice jako (m)
(m)
V (m) = diag[V i ], Λ(m) = diag[Λi ], D = diag[B i ] . Vektory q, f E a f I mají dimenzi n =
N X
ni a vektor x má dimenzi m =
i=1
N X
mi . Matice
i=1
Λ(m) je řádu m, matice D je řádu n a matice V (m) je obdélníková, typu (n, m). Pomocí zavedených vektorů a matic přepíšeme (40) do tvaru q = V (m) x 18
(43)
a (42) do tvaru ¨ + (V (m) )T DV (m) x˙ + Λ(m) x = (V (m) )T [f E + f I ] . Emx
(44)
Protože vazby mezi subsystémy jsou lineární diskrétní, existuje zřejmě vazební matice tuhosti K V a vazební matice tlumení B V tak, že lze potenciální energii kumulovanou v deformovaných vazbách vyjádřit ve tvaru kvadratické formy EpV = 21 q T K V q a ˙ Odtud Rayleighovu disipativní funkci lze analogicky vyjádřit ve tvaru RdV = 12 q˙ T B V q. plyne, že vektor silových účinků přenášených vazbami lze psát jako f I = −B V q˙ − K V q .
(45)
Dosazením (43) do (45) dostaneme f I = −B V V (m) x˙ − K V V (m) x . Dosazením tohoto výrazu do (44) dostaneme (m) T (m) ˙ ˙ ¨ +(V (m) )T DV (m) x+Λ ) K V V (m) x . x = (V (m) )T f E −(V (m) )T B V V (m) x−(V Emx
Sdružením členů s derivacemi stejného řádu vektoru hlavních modálních souřadnic získáme konečný tvar redukovaného modelu systému jako ¨ + (V (m) )T (D + B V )V (m) x˙ + [Λ(m) + (V (m) )T K V V (m) ]x = (V (m) )T f E . (46) E mx Redukovaný model v prostoru dimenze m má tedy známý tvar ¨ + B R x˙ + K R x = f E E mx R, kde redukovaná matice hmotnosti je jednotková, redukované matice tlumení a tuhosti mají tvar B R = (V (m) )T (D + B V )V (m) , K R = Λ(m) + (V (m) )T K V V (m) (m) T E a redukovaný vektor buzení je f E ) f . R = (V
Poznámka: Je-li původní model konzervativní (tedy bez tlumení a budících účinků), je i redukovaný model konzervativní tvaru ¨ + [Λ(m) + (V (m) )T K V V (m) ]x = 0 . Emx
(47)
Jeho vlastní frekvence aproximují vybrané vlastní frekvence původního modelu. Jestliže y i je (Eukleidovou normou normovaný) vlastní vektor redukovaného modelu příslušející jeho i−té vlastní frekvenci (tedy vektor ”malé” dimenze m), je zřejmě podle (43) vi = V (m) y i
(48)
zmíněné vlastní frekvenci příslušející vlastní vektor ”velké” dimenze n, jenž aproximuje vlastní vektor této frekvenci příslušející v původním modelu. Poznamenejme dále, že pokud potřebujeme znát pouze souřadnice i−tého vlastního vektoru, které popisují pouze (j) konfiguraci j−tého subsystému (označme tento subvektor vektoru v i jako v i ), stačí znát pouze subvektor 19
(j)
y i = [yJ+1,i , . . . , yJ+mj ,i ]T ; J =
j−1 X
ml
l=1
i−tého vlastního vektoru redukovaného modelu, přičemž podle (43) platí (j)
(m)
(j)
v i = V j y i ; j = 1, . . . , N . Této vlastnosti využíváme zejména v případě, kdy systém je složen ze subsystémů, z nichž jeden má pouze malý počet stupňů volnosti a řešitel potřebuje tvary kmitání právě tohoto subsystému. Dochází tím k velké úspoře paměti při výpočtu. Z hlediska úspory paměti při výpočtu i z hlediska výpočtového času je rovněž důležité, že k použití metody je nutné modální analýzu provádět pouze pro subsystémy a nikoliv pro celý systém. Příklad : Mějme soustavu tvaru přímého ”zleva” vetknutého řetězce s deseti stupni volnosti podle obr.10. Jeho parametry nechť jsou mi = i [kg] , i = 1, . . . , 10 ; kj = (11 − j) · 104 [N/m] , j = 1, . . . , 10 . Rozdělíme systém na dva subsystémy, kdy první subsystém je ”zleva” vetknutý řetězec se šesti stupni volnosti a druhý subsystém pak izolovaný řetězec se čtyřmi stupni volnosti. Provedeme redukce počtu stupňů volnosti systému uvažováním různého počtu hlavních tvarů kmitání subsystémů a budeme srovnávat získané vlastní frekvence (popřípadě i vlastní vektory) s příslušnými veličinami původního modelu. k1
m1
k2
m2
k3
m3
k4
k5
m4
m5
k6
m6
k7
m7
k8
1. subsystem
m8
k9
m9
k 10
m10
2. subsystem
Obrázek 9:
Řešení: Abychom mohli provádět výše popsaná srovnání, provedeme nejprve modální analýzu plného (”velkého”) modelu. Jeho vlastní frekvence jsou Ω1 = 12.903, Ω2 = 32.211, Ω3 = 52.546, Ω4 = 74.551, Ω5 = 99.238, Ω6 = 128.21, Ω7 = 164.35, Ω8 = 213.68, Ω9 = 292.26, Ω10 = 470.3, [rad/s] . M −normou normované vlastní vektory jsou
v1 =
0.011195 0.023614 0.037486 0.053073 0.070668 0.090607 0.11327 0.13908 0.16854 0.2022
; v2 =
0.024159 0.050772 0.079291 0.10841 0.13489 0.15267 0.15113 0.11249 0.0078444 −0.20878
; v3 =
20
−0.037363 −0.077731 −0.11778 −0.14961 −0.15921 −0.12677 −0.033719 0.11208 0.20699 −0.11754
; v4 =
0.05238 0.10735 0.15427 0.17114 0.12742 0.0041344 −0.15342 −0.16453 0.18458 −0.040497
;
v5 =
−0.071071 −0.14226 −0.18733 −0.15976 −0.022717 0.16411 0.15522 −0.21332 0.074193 −0.0083851
; v6 =
0.096521 0.18614 0.21046 0.089985 −0.14918 −0.19095 0.22768 −0.08746 0.014912 −0.0009659
; v7 =
0.13496 0.24441 0.20249 −0.079842 −0.26545 0.22887 −0.080582 0.014721 −0.001384 5.3205e − 5
;
−0.89472 −0.34832 0.2018 0.31001 −0.40476 0.32364 −0.048914 0.39606 0.091264 0.0045586 −0.13856 −0.3529 −0.00027519 0.026721 0.20317 v8 = −0.057235 ; v 9 = −0.0031853 ; v 10 = 1.1132e − 5 −3.0133e − 7 0.0092777 0.000242 5.2579e − 9 −1.1473e − 5 −0.00088573 −5.3635e − 11 3.1067e − 7 4.6528e − 5 2.4359e − 13 −3.6802e − 9 −1.0418e − 6 Nyní provedeme modální analýzu obou subsystémů. Vlastní frekvence prvního tému jsou
.
subsys-
(1)
(1)
(1)
(1)
(1)
(1)
Ω1 =31.724; Ω2 =98.541; Ω3 =158.04; Ω4 =213.41; Ω5 =292.26; Ω6 =470.3 [rad/s] a příslušnou maticí hmotnosti normované vlastní vektory pak
(1)
v1
=
−0.043661 −0.091686 −0.14341 −0.19633 −0.2449 −0.27854
; v (1) = 2
0.20361 0.32681 0.093325 (1) (1) v4 = ; v5 = −0.35567 0.20036 −0.044873 Vlastní frekvence druhého subsystému (2)
(2)
0.1101 0.22056 0.29128 0.25089 0.041349 −0.25025
−0.34832 −0.40477 0.39606 −0.13855 0.026672 −0.0028835 jsou (2)
; v (1) = 3
; v (1) = 6
−0.15054 −0.27603 −0.24486 0.052859 0.31218 −0.15632
;
−0.89472 0.31001 −0.048914 0.0045586 −0.00027518 0.000010773
.
(2)
Ω1 = 0, Ω2 = 32.419, Ω3 = 62.228, Ω4 = 99.728 [rad/s] a příslušnou maticí hmotnosti normované vlastní vektory pak
−0.20933 0.19129 −0.18175 −0.1715 0.27646 0.01845 −0.13718 −0.1715 (2) (2) (2) (2) . ; v 4 = ; v 3 = ; v2 = v 1 = −0.094688 −0.26939 −0.012657 −0.1715 0.010585 0.093789 0.24835 −0.1715 21
"
#
4 −4 Vazební matice tuhosti K V je řádu deset a obsahuje zaplněné pouze hnízdo ·104 −4 4 na místech šestého a sedmého řádku a sloupce. Zahrneme nejprve do redukovaného modelu vlastní hodnoty subsystémů odpovídající nejnižším vlastním frekvencím cca do 40[rad/s]. Bereme tedy m1 = 1 a m2 = 2. Redukovaný model bude mít tři stupně volnosti a bude mít tvar (47), kde
Λ
(m)
(1)2
Ω 1 = 0 0
0
0 0
(2)2 Ω1
(2)2 Ω2
0
; V
(m)
=
"
(1)
v 1 06,1 06,1 (2) (2) 04,1 v 1 v 2
#
.
Vlastní frekvence takto zkonstruovaného redukovaného modelu jsou Ω1 = 13.632 ; Ω2 = 32.213 ; Ω3 = 80.215[rad/s] . Srovnáním se třemi nejnižšími vlastními frekvencemi plného modelu vidíme, že aproximace je rozumně přesná pro první dvě vlastní frekvence (jež spadají do výše určeného rozsahu). Druhá vlastní frekvence je určena prakticky přesně a první s relativní chybou 5.8%. Třetí vlastní frekvence je (samozřejmě) nepoužitelná, neboť je určena s obrovskou relativní chybou (52.7%). Aproximace vlastních vektorů se kvantifikuje pomocí matice, vyjadřující tzv. MAC-kriterium. Jestliže v i je i−tý vlastní vektor plného modelu a ˆ i odpovídající vlastní vektor (dimenze n) vyčíslený pomocí redukovaného modelu, má v matice MAC prvky ˆs )2 ; r = 1, . . . , m , s = 1, . . . , m . M ACrs = (v Tr M v
(49)
MAC kriterium tedy ukazuje, nakolik vlastní vektory, určené prostřednictvím redukovaného modelu, tvoří, při postupném zastupování vlastních vektorů plného modelu, ortogonální normovaný systém vektorů prostoru Rn . Jestliže by byly ”dlouhé” vlastní vektory, určené prostřednictvím redukovaného modelu, určeny zcela přesně, byla by matice MAC jednotkovou maticí. Na obr.10 je zobrazena matice MAC pro výše zmíněný případ redukce. Je z tohoto obrázku vidět, že rovněž první dva vlastní vektory jsou určeny s rozumnou chybou, zatímco třetí vlastní vektor nikoliv.
1 0.8 0.6 0.4 0.2 0 1
2
3 2 3 1
Obrázek 10: 22
Zahrneme nyní do redukovaného modelu vlastní hodnoty subsystémů odpovídající nejnižším vlastním frekvencím cca do 100[rad/s]. Bereme tedy m1 = 2 a m2 = 4. Redukovaný model bude mít šest stupňů volnosti a bude mít tvar (47), kde
=
Λ(m)
V
(m)
(1)2
Ω1 0 0 0 0 0
=
"
0
0 0
(1)2
Ω2 0 0 0 0
0 0 0
(2)2
Ω1 0 0 0
0 0 0 0
(2)2
Ω2 0 0
(2)2
Ω3 0
(2)2
Ω4
(1)
(1)
0 0 0 0 0
v 1 v 2 06,1 06,1 06,1 06,1 (2) (2) (2) (2) 041 041 v 1 v 2 v 3 v 4
#
,
.
Vlastní frekvence takto redukovaného modelu jsou
Ω1 =, 12.936 Ω2 = 32.211, Ω3 = 52.697, Ω4 = 74.906, Ω5 = 99.239, Ω6 = 132.16 [rad/s] . Srovnáním se šesti nejnižšími vlastními frekvencemi plného modelu vidíme, že vypočítaná aproximace vlastních frekvencí je rozumná pro všech šest nejnižších vlastních frekvencí (včetně šesté vlastní frekvence, jež svou velikostí už vybočuje z předpokládaného frekvenčního intervalu). Druhá a pátá vlastní frekvence je aproximována prakticky přesně, ostatní s relativní chybou řádu desetin procenta a šestá vlastní frekvence je určena s relativní chybou tři procenta. Rovněž vlastní vektory příslušející k šesti nejnižším vlastním frekvencím jsou určeny s uspokojivou přesností, jak vyplývá ze znázornění MAC-kriteria pro tento případ uvedeného na obr.11. Rozdíl od jednotkové matice nelze vizuálně z obrázku vůbec vysledovat.
1 0.8 0.6 0.4 0.2 0 1 2 3
7 6
4 5
5
4
6
3 7
2 1
Obrázek 11: 23
Jestliže zahrneme do výpočtu ještě třetí vlastní frekvenci a odpovídající vlastní vektor prvního subsystému, bude aproximace odpovídat frekvenčnímu intervalu cca do 160 rad/s. Redukovaný model má pak sedm stupňů volnosti a vlastní frekvence jsou Ω1 =, 12.904 Ω2 = 32.211, Ω3 = 52.552, Ω4 = 74.565, Ω5 = 99.238, Ω6 = 128.29 Ω7 = 164.41 [rad/s] . Srovnáním se sedmi nejnižšími vlastními frekvencemi plného modelu je patrno, že tyto jsou uvedeným redukovaným modelem všechny aproximovány s minimálními relativními chybami nepřesahujícími desetinu procenta. Jestliže zahrneme do redukce vlastní hodnoty příslušející všem vlastním frekvencím obou subsystémů, je ”redukovaný” model přesně totožný s plným modelem. Na obrázku 12 je znázorněno MAC-kriterium shody vlastních vektorů obou modelů, kdy graf na tomto obrázku znázorňuje přesně jednotkovou matici desátého řádu.
1.4 1.2 1 0.8 0.6 0.4 0.2 0 1 2 3 4 5 6 7 8 9 10 1
2
3
4
5
6
7
8
9
10
Obrázek 12:
Poznámka: Z výrazu (47) je patrno, že změna modálních vlastností při zvedání počtu do redukce zahrnutých tvarů kmitání subsystémů je tím větší, čím větší (v absolutní hodnotě) jsou vazební souřadnice vlastních vektorů subsystémů. V našem případě se jedná o šestou souřadnici vlastních vektorů prvního subsystému a první souřadnici vlastních vektorů druhého subsystému (celkově tedy sedmou souřadnici). Zejména šesté souřadnice prvního subsystému jsou pro tři nejvyšší tvary kmitání tohoto subsystému řádově menší než než pro tři nejnižší tvary kmitání. Proto pro zmenšení chyb sedmi nejnižších vlastních frekvencí plného modelu nemá praktický význam zahrnování zmíněných vysokých frekvencí do redukovaného modelu. 24
6
Přibližné zahrnutí vedlejších tvarů kmitání subsystémů
Provedením modální analýzy přidruženého konzervativního modelu i−tého sybsystému získáme známé spektrální a modální (maticemi hmotnosti normované) matice Λi a V i , i = 1, . . . , N . Tvary kmitání každého subsystému rozdělíme (nejčastěji podle velikosti příslušných vlastních frekvencí) na mi hlavních tvarů, si vedlejších tvarů a ostatní tvary. Hlavní tvary ovlivní redukovaný model zásadním způsobem podobně jako v minulé kapitole, vedlejší tvary jej ovlivní korekčně (nepřidají už počet stupňů volnosti redukovaného modelu) a ostatní tvary neovlivní redukovaný model vůbec. Matice Λi a V i (po eventuální permutaci pořadí vlastních tvarů kmitání a k nim příslušných vlastních frekvencí) lze psát ve tvaru rozdělených matic (m)
(s)
(o)
(m)
(s)
(o)
Λi = diag[Λi , Λi , Λi ] ; V i = [V i , V i , V i ] , i = 1, . . . , N . (m)
(s)
(50) (o)
Poznamenejme, že matice Λi je řádu mi , matice Λi je řádu si a matice Λi je řádu (m) (s) ni − mi − si . Matice V i je obdélníková typu (ni , mi ), matice V i je obdélníková typu (o) (ni , si ) a matice V i je typu (ni , ni − mi − si ). Provedeme nyní pro každý subsystém modální transformaci (m)
(m)
q i = V i xi (m)
(s)
(s)
+ V i xi ,
(51)
(s)
kde xi jsou hlavní a xi vedlejší modální souřadnice i−tého subsystému. Dosazením (51) do modelu (38) vznikne (m) (m) (m) (m) ˙ (s) ˙ i + B i V (s) ¨ (s) ¨ i + M i V (s) M iV i x i + i x i + BiV i x i x (m)
(m)
K i V i xi
(s)
(s)
I + K i V i xi = f E i + fi . (m)T
(s)T
Násobíme tuto rovnici zleva nejprve maticí V i a poté i maticí V i hlavních a vedlejších tvarů jsou disjunktní, platí zřejmě (m)T
Vi
(m)T
Vi
(s)
(s)T
M iV i
(s)
(s)T
K iV i
M i V i = 0mi ,si ; V i K i V i = 0mi ,si ; V i
(m)
(m)
. Protože množiny
= 0si ,mi ; = 0si ,mi ;
a tudíž z předchozí rovnice obdržíme (m)
(m) (m) (m)T (m)T (s) (s) (m) (m) I B i V i x˙ i +V i B i V i x˙ i +Λi xi =V i [f E i +f i ] ,
(52)
(s) (s) (s)T (s)T (s) (s) (m) (m) I B i V i x˙ i + V i B i V i x˙ i + Λi xi = V i [f E i + fi ] .
(53)
(m)T
¨ i +V i E mi x (s)
(s)T
¨i + V i E si x
(s)
(s)
(m)T
V rovnici (53) zanedbáme členy s derivacemi vektoru xi a člen V i B i V i (m)T (s) (a tím je zanedbán i v rovnici (52) člen V i B i V i ) Rovnice (52) a (53) mají potom jednodušší tvar (m)
¨i E mi x
(m)T
+Vi
(m) (m) (m)T (m) (m) I B i V i x˙ i + Λi xi = V i [f E i = fi ] , (s)
(s)
(s)T
Λi xi = V i
I [f E i + fi ] .
(54) (55)
Za předpokladu, že žádný vedlejší tvar i−tého subsystému nemá nulové vlastní číslo, vyplývá z (55) 25
(s)
(s)−1
xi = Λi
(s)T
Vi
I [f E i + fi ] .
Dosazením tohoto výrazu do modální transformace (51) máme (m)
(m)
q i = V i xi
(s)
(s)−1
+ V i Λi
(s)T
Vi
I [f E i + fi ] .
(56)
Zaveďme nyní globální vektory celého systému jako (m)T
q T = [q T1 , . . . , q TN ], xT = [x1
(m)T
, . . . , xN
],
E T T I I T I T T T f E = [(f E 1 ) , . . . , (f N ) ] , f = [(f 1 ) , . . . , (f N ) ]
a matice jako (m)
(m)
V (m) = diag[V i ], Λ(m) = diag[Λi ], D = diag[B i ] , (s)
(s)
V (s) = diag[V i ], Λ(s) = diag[Λi ] . Vektory q, f
E
I
a f mají dimenzi n =
N X
ni a vektor x má dimenzi m =
mi . Matice
i=1
i=1
Λ(m) je řádu m, matice Λ(s) je řádu s =
N X
N X
si , matice D je řádu n, matice V (m) je
i=1
obdélníková typu (n, m) a matice V (s) je typu (n, s). Pomocí zavedených vektorů a matic přepíšeme (54) a (56) do tvaru ¨ + V (m)T DV (m) x˙ + Λ(m) x = V (m)T [f E + f I ] , Emx
(57)
q = V (m) x + V (s) Λ(s)−1 V (s)T [f E + f I ] .
(58)
Podobně jako v předchozí kapitole, existují vazební matice tlumení B V a tuhosti K V (obě ”velkého” řádu n), že pro vazební síly platí f I = −B V q˙ − K V q .
(59)
Zanedbáním vlivu vedlejších tvarů kmitání na zobecněné rychlosti lze derivací ˙ Dosazením (58) za q a posledního výrazu za q˙ do (59) dostaneme (58) psát q˙ = V (m) x. f I = −B V V (m) x˙ − K V V (m) x − K V V (s) Λ(s)−1 V (s)T [f E + f I ] . Označíme-li G = V (s) Λ(s)−1 V (s)T ,
(60)
což je matice ”velkého” řádu n, dostaneme odtud f I = −(E n + K V G)−1 [B V V (m) x˙ + K V V (m) x + K V Gf E ] . Dosazením tohoto vztahu do (57) vznikne ¨ + V (m)T DV (m) x˙ + Λ(m) x = V (m)T f E − Emx −V (m)T (E n + K V G)−1 [K V V (m) x + B V V (m) x˙ + K V Gf E ] . Sdružením členů se stejnými derivacemi vektoru hlavních modálních souřadnic dostáváme redukovaný model ve tvaru 26
¨ + [V (m)T DV (m) + V (m)T (E n + K V G)−1 B V V (m) ]x+ ˙ Emx +[Λ(m) + V (m)T (E n + K V G)−1 K V V (m) ]x = V (m)T [E n − (E n + K V G)−1 K V G]f E . Redukovaný model v prostoru dimenze m má tedy známý tvar ¨ + B R x˙ + K R x = f E E mx R,
(61)
kde redukovaná matice hmotnosti je jednotková, redukované matice tlumení a tuhosti mají tvar B R = V (m)T [D + (E n + K V G)−1 B V ]V (m) , K R = Λ(m) + V (m)T (E n + K V G)−1 K V V (m)
(62)
(m)T a redukovaný vektor buzení je f E [E n − (E n + K V G)−1 K V G]f E . R = V Redukovaný model je zde stejné dimenze jako v předchozí kapitole. Korekční příspěvek vedlejších tvarů kmitání subsystémů je zaveden pomocí matice G v (60). Jestliže se žádný vedlejší tvar neuvažuje, je matice G nulovou maticí řádu n. Proto v tomto případě (E N + K V G)−1 = E n a model (61) přejde na analogický model z předchozí kapitoly.
Poznámka: Je-li původní model konzervativní (tedy bez tlumení a budících účinků), je i zde redukovaný model konzervativní tvaru ¨ + [Λ(m) + (V (m) )T (E n + K V G)−1 K V V (m) ]x = 0 . E mx
(63)
Jeho vlastní frekvence aproximují vybrané vlastní frekvence původního modelu. Jestliže y i je (Eukleidovou normou normovaný) vlastní vektor redukovaného modelu příslušející jeho i−té vlastní frekvenci (tedy vektor ”malé” dimenze m), transformuje se na vlastní vektor v i ”velké” dimenze podle (58), kam za f I dosazujeme (59) při uvažování konzervativní soustavy (tedy pro B V = 0 a f E = 0). Odtud dostáváme transformaci v i = V (m) y i + V (s) Λ(s)−1 V (s)T (−K V v i ) . Dosazením do tohoto výrazu matici G podle (60), dostáváme odtud po drobné úpravě vztah mezi ”krátkými” vlastními vektory y i a jejich ”dlouhými” protějšky v i ve tvaru v i = (E n + GK V )−1 V (m) y i .
(64)
Tyto vektory aproximují vlastní vektory původního modelu. Poznamenejme ještě, že zatímco v koeficientových maticích redukovaného modelu (62) jest nutno invertovat (velkou) matici E n + K V G, v přepočtu vlastních vektorů (64) se jedná o matici E n + + GK V . Poznámka: Výhoda korekce modelu vedlejšími tvary kmitání je na první pohled eliminována nutností invertovat dvě velké matice řádu rovného dimenzi původního modelu. Naštěstí však vazební matice tuhosti K V je zpravidla velmi řídká. Jestliže například systém se skládá ze dvou subsystémů, jež jsou spojeny jedinou vazbou, má zmíněná matice nenulové prvky pouze ve hnízdě řádu dva, propojujícím souřadnice plného modelu, jež vazba spojuje. Lze tedy s výhodou využít Hauseholderovy identity (viz příslušné téma) a invertovat pouze matici takového řádu, kolik činí počet nenulových řádků (sloupců) vazební matice K V . V předchozím příkladě by tedy stačila pro výpočet inverze matice řádu dva. Kromě toho je možno uspořit paměť počítače, protože si nemusíme pamatovat celou matici G-viz rovněž téma ”Hauseholderova identita”. 27
