Analiza modalna z wykorzystaniem metody elementów skończonych
Analizowanym obiektem jest konstrukcja wsporcza morskiej turbiny wiatrowej (MTW) w skali laboratoryjnej (1:40). W pierwszym etapie przeprowadzono numeryczną analizę modalną konstrukcji wsporcze MTW, a następnie przeprowadzono korelację uzyskanych wyników w zakresie wielkości wartości własnych, jak i postaci drgań struktury. Wyniki uzyskano dla warunków podparcia w konfiguracji swobodnego brzegu. W drugim etapie celem była ocena skorelowanych wyników uzyskanych na drodze eksperymentalnej analizy modalnej z wynikami symulacji komputerowej dla struktury podpartej na elementach podatnych, oraz określenie jakości uzyskiwanych wyników symulacji komputerowej w warunkach niepewnych parametrów podparcia struktury. W przypadku obiektu rzeczywistego podparcie zrealizowano stosując gumowe pasy, oraz gumowe bloki na których ustawiono konstrukcję trójnoga MTW. Ostatnim działaniem była próba oszacowania podatności posadowienia MTW w przypadku wykorzystania gumowych bloków.
Pierwotny model MES
Pierwszy model MES powstał przed fizyczną realizacją konstrukcji modelu testowego w oparciu o geometrię bryłową dostarczoną przez CTO S.A. (Rysunek 20).
Rysunek 20. Dostarczony przez CTO S.A. model bryłowy trójnoga.
Ze względu na nieznaczną grubość ścianek trójnoga zdecydowano się na modelowanie z wykorzystaniem elementów płytowych. Model bryłowy zastąpiono modelem powierzchniowym. Dodatkowo wykorzystano symetrię modelu (Rysunek 21).
Rysunek 21. Geometria powłokowa trójnoga powstała na podstawie modelu bryłowego.
Siatkę MES modelu zaktualizowanego zaprezentowano na (Rysunku 22). Siatka MES składała się z: 48 elementów trójkątnych, 7108 elementów czworobocznych i 7027 węzłów. Obiekt rzeczywisty wykonany został z aluminium, w celu budowy modelu MES założono następujące parametry: moduł Young’a E=70GPa, liczba Poisson’a ratio n=0.33, gęstość masy r=2750kg/m3. W pierwszym modelu przypisano identyczny parametr własności elementu skończonego dla powtarzalnych składowych części konstrukcji.
Rysunek 22. Siatka pierwotnego modelu MES.
Wyniki analizy modalnej pierwotnego modelu MES
Postacie i częstości drgań własnych przedstawiono w Tabeli 7. Ze względu na założoną symetrię pojawiają się one w parach (dwie identyczne wartości głównych centralnych masowych momentów bezładności I1=6.015kgm2, I2=I3=13.806kgm2 – osie oznaczone zgodnie z Rysunkiem 23). Jądrem obliczeniowym był MSC/Nastran. Pre i postprocesorem LMS Virtual.Lab.
Tabela 7. Wyznaczone postacie drgań własnych i odpowiadające im częstości pierwotnego modelu MES.
w1 = 62.9Hz |
w2 = 62.9Hz |
w3 = 104.1Hz |
w4 = 104.1Hz |
w5 = 221.2Hz |
w6 = 252.6Hz |
w7 = 252.6Hz |
w8 = 262.1Hz |
w9 = 284.2Hz |
w10 = 284.2Hz |
Rysunek 23. Układ osi głównych centralnych.
Korelacja modelu MES z modelem pomiarowym.
Po uzyskaniu wyników testu obiektu rzeczywistego mierzonego w warunkach swobodnego brzegu, oba modele poddano porównaniu. Kryterium zgodności stanowił MAC (Modal Assurance Criterion – modalne kryterium zgodności), obliczane zgodnie ze wzorem 3:
(3)
gdzie:
Ψitest – wektor postaci własnej uzyskanej na drodze pomiaru
ΨjFE – wektor postaci własnej uzyskanej w trakcie symulacji.
Wartość MAC może przyjmować wartości z zakresu 0 do 1, 1 – postacie są identyczne, 0 – zupełny brak korelacji w postaci drgań własnych. Ze względu na różnice w ilości węzłów w modelu MES i modelu estymowanym na drodze testu (tylko 55 węzłów), model MES musiał zostać zredukowany. Do porównania brano węzły modelu MES leżące najbliżej węzłów modelu testowego. Prezentują to żółte punkty na Rysunku 24. Dodatkowo, ze względu na braku rotacyjnych stopni swobody w węzłach modelu testowego, porównywano jedynie translacyjne stopnie swobody.
Rysunek 24. Nałożone na siebie geometrie modelu MES i modelu pomiarowego.
Wartości wielkości MAC przedstawiono na rysunku 25. Dodatkowo w Tabeli 8 podano wartości częstości drgań własnych wraz z błędem względnym.
Rysunek 25. Macierz MAC porównania wyników eksperyment/symulacja.
Tabela 8. Porównanie wyniku częstości uzyskanej podczas symulacji z wynikami eksperymentalnej analizy modalnej.
| FE ω | EMA ω | error [%] | |
| 1 | 62.9 | 66.13 | 4.8843188 |
| 2 | 62.9 | 69.24 | 9.1565569 |
| 3 | 104.1 | 91.12 | 14.244952 |
| 4 | 104.1 | 94.86 | 9.7406705 |
| 5 | 221.2 | 210.12 | 5.2731772 |
| 6 | 252.6 | 247.92 | 1.8877057 |
| 7 | 252.6 | 258 | 2.0930233 |
| 8 | 262.1 | 267.1 | 1.8719581 |
| 9 | 284.2 | 283 | 0.4240283 |
| 10 | 284.2 | 287.4 | 1.1134308 |
Problemy występują w parach 3-3 i 8-8. Może być to spowodowane nieuwzględnieniem naprężenia w modelu numerycznym naprężenia od ciężaru własnego konstrukcji, który powoduje podwyższenie wartości sztywności. Dodatkowym elementem może być umieszczeni wzbudników bezpośrednio na gruncie, co powoduje znaczny wpływ popychaczy w trakcie testu na jego wyniki. Konieczne okazuje się też przypisanie oddzielnych własności elementów skończonych dla powtarzalnych części konstrukcji wsporczej (tj. rezygnacja z założonej symetrii).
Analiza wrażliwości
W celu poprawy modelu numerycznego stosowana jest procedura strojenia z wykorzystaniem metod optymalizacyjnych. Aby wybrać jak najkorzystniejszy zbiór parametrów modelu do optymalizacji, przeprowadza się analizę wrażliwości. Uzyskane informacje podają jak ilościowa zmiana parametru modelu wpływa na jakościową poprawę korelacji modelu MES i eksperymentalnego. Określenie wpływu zmian można obliczyć zgodnie ze wzorami 4 i 5:
(2) 
(3)
gdzie:
[M] – macierz mas,
[K] – macierz sztywności,
[F] – wektor postaci drgań własnych,
m – masa modalna,
w – wartość częstości drgania własnego,
q – testowany parametr modelu.
W trakcie analizy okazało się, że nie będzie możliwe dostrojenie modelu w którym założono symetrię parametrów modelu MES (Rysunek 26). Przyjęty określający grubość nogi trójnoga uniemożliwia oddzielne strojenie częstości 8, 9 i 10. Sytuacja poprawia się gdy parametry grubości nóg zostały odseparowane (Rysunek 27).
Rysunek 26. Wyniki analizy wrażliwości dla wspólnych parametrów dla elementów powtarzalnych w konstrukcji.
Rysunek 27. Wyniki analizy wrażliwości dla odseparowanych parametrów dla elementów powtarzalnych w konstrukcji.
Zmiany w modelu MES
Ze względu na modyfikację laboratoryjnego modelu konstrukcji wsporczej morskiej turbiny wiatrowej zaistniała konieczność zaktualizowania modelu numerycznego (MES) struktury. Modyfikacja obiektu rzeczywistego polegała na dodaniu otworów zalewowych, oraz kołnierza mającego służyć symulacji procesu pękania na jednym z górnych tężników konstrukcji (Rysunek 28).
Rysunek 28. Modyfikacja konstrukcji laboratoryjnego modelu wsporczej morskiej turbiny wiatrowej. Dodane otwory zalewowe, oraz kołnierz na jednym z górnych tężników.
Zmieniony został model geometryczny (Rysunek 29). Siatkę MES modelu zaktualizowanego zaprezentowano na (Rysunek 30). Pierwotna siatka MES składała się z: 48 elementów trójkątnych, 7108 elementów czworobocznych i 7027 węzłów. Po modyfikacji siatka składa się z: 390 elementów trójkątnych, 32470 elementów czworobocznych i 32573 węzłów. Masa struktury pierwotnej wynosi 32.011kg, masa struktury wsporczej po modyfikacji wynosi 32.499kg. Podobnie jak w modelu pierwotnym każdy z elementów strukturalnych ma oddzielne własności materiałowo-geometryczne (moduł Young’a, liczbę Poissona, grubość). Umożliwia to strojenie każdego z elementów.
Rysunek 29. Geometria powłokowa trójnoga z dodanymi otworami zalewowymi i kołnierzem.
Rysunek 30. Topologia siatki MES po wprowadzeniu zmian.
Podobnie jak w przypadku modelu pierwotnego MES przeprowadzono korelację z modelem eksperymentalnym, analizę wrażliwości, oraz strojenie. W wyniki przeprowadzonych działań uzyskano wyniki przedstawione na Rysunkach 31 i 32.
Rysunek 31. Porównanie zgodności modelu MES i modelu testowego.
Rysunek 32. Względna różnica częstości drgań własnych modelu MES i modelu testowego.
Uzyskane wyniki wskazują na dobre odwzorowanie dynamiki struktury poprzez model MES w porównaniu do modelu testowego, zarówno w odniesieniu do postaci, jak i częstości drgań własnych obiektu dla konfiguracji swobodnego brzegu.
Próby korelacji i dostrajania modelu MES z testem wykonanym dla posadowienia podatnego
Analizując wyniki korelacji wnioskowano, że podatny sposób podparcia powoduje znaczące problemy ze skojarzeniem wyników analizy numerycznej z wynikami eksperymentalnymi. Podczas eksperymentu warunki podatnego podparcia symulowano poprzez posadowienie konstrukcji trójnoga na blokach, oraz taśmach wykonanych z gumy. Szczególnie trudnym zadaniem było skojarzenie wyników symulacji i testu w przypadku posadowienia na gumowych pasach. Może to być spowodowane łatwością zmiany geometrii pasa gumy, a co za tym idzie sztywności podparcia, ze względu na interakcję gumowego pasa w trakcie wymuszania struktury, jak i możliwością przemieszczenia względnego samostykającej się taśmy w trakcie przygotowywania kolejnych pomiarów. Z powyższych względów doradza się nie stosowanie tego typu podparcia i nie brania go pod uwagę w symulacji korelacji modeli.
W przypadku podparcia struktury wsporczej na blokach gumowych, porównanie korelacji modeli w warunkach swobodnego brzegu i skrępowanej, również nie przynosiło zadowalającego efektu. Problemem było ustalenie sztywności wzdłużnej gumowych bloków w trzech kierunkach pionowych i sześciu kierunkach poziomych (bloki miały nieznacznie różniącą się geometrię). Dopiero przybliżenie systemu modelem oscylatora harmonicznego umożliwiło nieznaczną poprawę sytuacji (Rysunek 33). Trójnóg traktowany był jako masa skupiona (32,5kg), natomiast bloki gumowe zastąpiono elementem sprężystym. Wartości podatności szacowano dla częstości własnej 120Hz, uzyskanej w trakcie testu
Rysunek 33. Uzyskana najlepsza korelacja i różnica częstości przy wykorzystaniu oszacowanej sztywności pionowej podatnego podparcia (założono, że bloki są ortotropowe, ale identyczne).
W celu sprawdzenia możliwości przyszłego stosowania podparcia podatnego w badaniach i analizie korelacji modeli, numerycznego i testowego, dokonano modyfikacji modelu MES. Założono, że sztywność wzdłużna na każdym kierunku każdego z bloków jest różna (9 dodatkowych zmiennych projektowych). Dodatkowo założono, że aluminiowy trójnóg można uznać za bryłę sztywną w stosunku do podatności bloków gumowych. Umożliwiło to nieuwzględnianie sztywności rotacyjnych trójwymiarowego elementu podatnego o sześciokrotnej ogólnej liczbie swobody (element klasy Bushing w nomenklaturze NASTRAN). Numeracje elementów podatnych łączących trójnóg z nieodkształcalnym podłożem zaprezentowano na (Rysunek 34).
Rysunek 34. Numeracja elementów podatnych.
B1 B2 B3
Do oszacowania wpływu sztywności elementów podatnych na zmianę wartości częstości własnej modelu numerycznego (MES), wykorzystano metodę DOE (Design of Experiment). Metoda umożliwia śledzenie wpływu zmian wartości wejściowych na zmianę wartości wyjściowych, poprzez wprowadzenie dyskretnego zbioru wartości parametru wejściowego, przy czym zmieniany może być więcej niż jeden parametr. Obliczenia dokonano na KDM (Komputerach Dużej Mocy) w TASK (Trójmiejskiej Akademickiej Sieci komputerowej) z wykorzystaniem oprogramowania Noesis Optimus 10.13 i MSC. Nastran.
Ze względu na znaczący wpływ sztywności pionowej elementów podatnych, zakres zmian numerycznych wartości tej wielkości został ustalony w wąskim przedziale, przeciwnie do zakresu zmian wartości sztywności poprzecznych.
Proces modelowania polega DOE polegał na określeniu: wsadowego pliku źródłowego (bushing_DOE.dat), źródłowego pliku wynikowego (bushing_doe.f06), wektora wielkości wejściowych – wraz ze wskazaniem miejsca występowania wymaganych parametrów w pliku (InputArray1 – 9 sztywności), wektora wielkości wynikowych – wraz ze wskazaniem miejsca występowania wymaganych parametrów w pliku (Frequen – 15 częstości), akcji umożliwiającej wygenerowanie pliku wyjściowego z pliku wejściowego (komendy uruchomieniowej programu MSC.Nastran), oraz wzajemnych relacji powyższych obiektów (Rysunek 35).
Rysunek 35. Schemat wykorzystania danych w programie Optimus.
Wartości nominalne parametrów, ich zakres zmienności, oraz zastosowane nazewnictwo parametrów wejścia/wyjścia zaprezentowano na (Rysunek 36), oraz (Rysunek 37). Sztywności oznaczono, odpowiednio: BK1,BK4, BK5 – sztywności w kierunku x (zgodnie z UW modelu) elementu podatnego (EP): 1, 2, 3; indeksy: 2, 5, 8 – oznaczają sztywność w kierunku osi y EP; natomiast 3, 6, 9 – oznaczają sztywności w kierunku osi pionowej EP.
Rysunek 36. Definicja wektora wielkości wejściowych.
Rysunek 37. Definicja wektora wielkości wejściowych.
W zaproponowanym podejściu wykorzystano model drugiego rzędu wykorzystujący pełnoczynnikową metodę DOE trzeciego poziomu (3 Level Full Factorial method). Ze względu na ilość wszystkich możliwych kombinacji parametrów modelu numerycznego, sumaryczna liczba eksperymentów wyniosła 19683 (Rysunek 38).
Rysunek 38. Definicja wykorzystanej metody 3 Level Full Factorial.
Ze względu na dużą ilość danych, oraz zależność analizowanego stanu od kilku parametrów (model drugiego rzędu), wygodnym narzędziem analizy rezultatów eksperymentów numerycznych jest metoda powierzchni odpowiedzi (RSM – Response Surface Model). Umożliwia ona trójwymiarową prezentacje wyników eksperymentu numerycznego i pokazuje w postaci powierzchni zmianę odpowiedzi systemu przy ustalonych dwóch wartościach parametru. (Rysunek 39), (Rysunek 40) i (Rysunek 41) przedstawiają przykładową zmianę w trzeciej częstości drgań własnych (odpowiadającą drganiom w osi pionowej obiektu – wyznaczonym z wykorzystaniem modelu oscylatora harmonicznego) dla ustalonych wartości sztywności EP 1 i 2 w kierunku pionowym, oraz zmianie pozostałych parametrów. Wyniki analizy wskazują na niewielką czułość częstości drgań własnych trójnoga na zmiany sztywności w kierunku poprzecznym do osi pionowej obiektu, jedyny realny wpływ ma pozostała trzecia sztywność EP 3 w kierunku osi pionowej. Obniżanie wartości sztywności poziomej EP powoduje, że trójnóg zaczyna drgać jak bryła sztywna. Natomiast zwiększanie sztywności w tych kierunkach powoduje pojawienie się częstości niewystępujących w wynikach pomiaru eksperymentalnego (zarówno co do wartości, jak i wektora własnego). Może to świadczyć o znikomym ograniczeniu możliwości ruchu trójnoga w kierunkach zawartych w płaszczyźnie jego posadowienia. Dodatkowo analizując wykresy rozproszenia, przedstawiające korelacje pomiędzy wartościami wyliczonymi w DOE, a wartościami estymowanymi przez model, widać, że zakres niezgodności jest zdecydowanie za duży. Należałoby znacznie ograniczyć zakres zmienności parametru, co skutkuje jeszcze mniejszym wpływem sztywności w kierunkach x i y modelu, na zmianę częstości drgań układu. Na (Rysunek 42) zaprezentowano przykładowy wykres rozproszenia dla przypadku częstotliwości 3. Sztywność w kierunku osi pionowej jest zdecydowanie za mała, różni się aż i 12 razy od wyznaczonej eksperymentalnie, co świadczy o tym, że sztywności gumowych bloków nie można wyznaczyć na podstawie znajomości wyłącznie: masy układu i częstości rezonansowej (pomimo, że postać przedstawia niezdeformowany trójnóg).
Rysunek 39. Zależność wartości 3 częstości własnej dla ustalonych wartościach sztywności pionowych: BK3 i BK6 (pozostałe sztywności o wartości minimalnej).
Rysunek 40. Zależność wartości 3 częstości własnej dla ustalonych wartości sztywności pionowych: BK3 i BK6 (pozostałe sztywności o wartości maksymalnej).
Rysunek 41. Zależność wartości 3 częstości własnej dla maksymalnych wartości sztywnościach pionowych.
Rysunek 42. Zależność wartości 3 częstości własnej dla maksymalnych wartości sztywnościach pionowych.
Wnioski
Modele: testowy i numeryczny, które wstępnie wykazywały dość dobrą korelację, stały się bardzo trudne do skorelowania, a co się z tym wiąże do zidentyfikowania, po zamianie warunków brzegowych z konfiguracji wolnego brzegu do konfiguracji sprężystego podparcia. Pomimo rezygnacji ze stosowania standardowych warunków brzegowych, tj. odbierania stopni swobody węzłom, na rzecz zastosowania podparcia podatnego w postaci elementów podatnych, w dalszym ciągu struktury wykazywały bardzo małe podobieństwo. Za taki stan odpowiada nieznajomość/niemożliwość wyznaczenia parametrów opisujących ilościowo własności strukturalne podpór. Pomimo prób dostrojenia elementów podatnych z wykorzystaniem jednej z postaci drgań własnych, bardzo przypominającej zachowanie oscylatora harmonicznego, nie udało się prawidłowo oszacować sztywności elementów podatnych w kierunku pionowym trójnoga. Wynikowa wartość częstości drgań trzeciej postaci własnej uzyskana w trakcie eksperymentów numerycznych oscylowała w pobliżu wartości 10Hz (przy wymaganych 120Hz). Dzięki przeprowadzonej analizie statystycznej, z wykorzystaniem DOE i RSM, możliwe było zauważenie trendów w zmianie odpowiedzi strukturalnej modelu wirtualnego, spowodowanych zmianą parametrów wejściowych. Na tej podstawie stwierdzono pomijalnie mały wpływ na trzecią wartość i postać własną innych sztywności, niż tylko te w kierunku pionowym. Dodatkowo DOE wskazało Kierunek poszukiwania bliższej fizyce zjawiska wartości parametrów podparcia. Wartości sztywności należy szukać znacznie powyżej oszacowanej wartości (wartości nominalnej) sztywności pionowej bloków (19700N/m) i znacznie poniżej maksymalnej wartości sztywności poprzecznych (1000N/m). Z analizy można również dowiedzieć się, czy zakres przeszukiwania pożądanych wartości nie jest zbyt rozległy i nie obarcza przez to analizy zbyt dużym uchybem.
Powyższe spostrzeżenia mogą stanowić ważną wskazówkę dla, np. technik SHM. W przypadku rzeczywistych morskich turbin wiatrowych, montowanych z posadowieniem na dnie, bardzo istotne będzie rozpoznanie własności strukturalnych dno-podpora, szczególnie w przypadku dużego gradientu w wartościach własności geologicznych na obszarze zajmowanym przez farmę wiatrową.
Przedstawiony przykład był znacznie mniej skomplikowany w analizie ze względu na otrzymane wcześniej zgodności modeli: testowego (w warunkach free-free) i numerycznego (modelowanego MES). Pojawiająca się niezgodność od razu została zdiagnozowana jako wpływ warunków podparcia, co ograniczyło obszar przeszukiwania jedynie do 9 parametrów. W przypadku braku pewności, pozostałe parametry konstrukcyjne związane wyłącznie z ze strukturą samego trójnoga, musiałyby być dodatkowo włączone do zbioru parametrów przeszukiwania. Taki stan rzeczy może poważnie utrudnić przeprowadzenie, lub nawet uniemożliwić, przeprowadzenie analizy statystycznej.
Literatura:
[1] Heylen, W., Lammens, S., Sas, P., Modal analysis theory and testing, Katholieke Universiteit Leuven, Belgia, 2009.
[2] Maciej Kahsin, Analiza i opracowanie korelacji wraz z oceną zgodności modeli Metody Elementów Skończonych modeli trójpodporowej konstrukcji wsporczej morskiej turbiny wiatrowej z pomiarami modelu konstrukcji w skali laboratoryjnej dla konfiguracji podatnego podparcia. Analiza wpływu zmienności pomiarów parametrów dynamicznych konstrukcji na dokładność dostrajania modeli numerycznych w ujęciu eksperymentalnym i teoretycznym, opracowanie pracy zleconej w ramach projektu AQUILO, PBS1/A6/8/2012, Gdańsk 2015.
[3] Maciej Kahsin, Marcin Łuczak, Numerical Model Quality Assessment of Offshore Wind Turbine Supporting Structure Based on Experimental Data, IWSHM2015, Stanford, USA, 2015.
[4] Maciej Kahsin, Marcin Łuczak, Bart Peeters, Use and assessment of preliminary FE model results within testing process of offshore wind turbine supporting structure, EURODYN 2014, Porto, Portugal 2014.
[5] Maciej Kahsin, Rozbudowa numerycznego modelu trójpodporowej konstrukcji wsporczej morskiej turbiny wiatrowej oraz analiza korelacji z pomiarami modelu konstrukcji w skali laboratoryjnej dla konfiguracji swobodnego brzegu, opracowanie pracy zleconej w ramach projektu AQUILO, PBS1/A6/8/2012, Gdańsk 2014.
[6] Dokumentacja techniczna, Optimus Theoretical Background, Leuven, Belgium 2014.
[7] Warren, C., Niezrecki, C., Avitabile, P., Comparison of FRF Measurements and Mode Shapes Determined using Optically Image Based, Laser, and Accelerometer Measurements, Mechanical Systems and Signal Processing, Volume 25, Issue 6, August 2011, Pages 2191–2202
[8] Bi, S., Ren, J., Wang, W., Elimination of Transducer Mass Loading Effects in Shaker Modal Testing, Mechanical Systems and Signal Processing, Volume 38, Issue 2, 20 July 2013, Pages 265–275
[9] Cakar, O., and Sanliturk, K., Elimination of Transducer Mass Loading Effects from Frequency Response Functions, Mechanical Systems and Signal Processing, Volume 38, Issue 2, 20 July 2013, Pages 265–275
[10] LMS International, The LMS Theory and Background Book—Analysis and Design, Manual of Test.Lab Revision 5, LMS International, (2005),
[11] Allemang, Randall J. The modal assurance criterion–twenty years of use and abuse, Sound and Vibration 37.8 (2003): 14-23.
[12] LMS International, The LMS Theory and Background Book—Analysis and Design, Manual of Virtual.Lab Revision 11, LMS International, (2011),
[13] Fox, R., L., Kapoor, M,. P., Rates of Change of Eigenvectors and Eigenvalues, AIAA Journal, Vol. 6, 1968, Pages 2426-2429.