Metoda Elementów Skończonych w Modelowaniu Układów Mechatronicznych

(1)

_________________________________________________________________

Metoda Elementów Skończonych w

Modelowaniu Układów Mechatronicznych

_________________________________________________________________

Układy belkowe (Scilab)

(2)

str.1

I. Podstawowe informacje

Dla elementów belkowych oprócz rozciągania (ściskania) występuje również zginanie oraz ścinanie.

W węzłach elementów belkowych wystąpią dwa stopnie swobody ^* (dla elementów prętowych był to tylko jeden stopień swobody – przemieszczenie wzdłuż osi x).

W każdym z węzłów wyróżniamy stopnie swobody przemieszczeniowe (indeksy nieparzyste) oraz przemieszczenia kątowe (indeksy parzyste).

Macierz sztywności płaskiego elementu belkowego ma postać:

Dla wektora obciążeń F indeksy nieparzyste odpowiadają siłom węzłowym, natomiast nieparzyste momentom skupionym w odpowiednich węzłach.

Agregacja, czyli łączenie lokalnych macierzy sztywności w globalny układ równań uwzględnia wszystkie stopnie swobody w węzłach.

W metodzie kary w zależności od typu podpory modyfikujemy odpowiednie elementy macierzy uwzględniające zarówno przemieszczenia jak i obroty w węzłach.

Po rozwiązaniu macierzowego równania MES

otrzymujemy wektor węzłowych przemieszczeń Q

Pozostałe wielkości, jak reakcje podporowe, odkształcenia oraz naprężenia w elementach wyznaczane są z odpowiednich formuł. W ramach tego laboratorium obliczenia MES dla belek dotyczyć będą wyłącznie wyznaczania przemieszczeń.

________________________________________________________________________________

(*) tak jest tylko w przypadku rozważanych płaskich elementów belkowych. W ogólnym przypadku przestrzennym, w którym tego typu elementami można modelować również ramy, w każdym z węzłów wystąpi 6 stopni swobody (3 przemieszczenia oraz 3 obroty)

3 2 3 2

2

3 2

12 6 12 6

4 6 2

12 6

4

e z

L L L L

L L L

EI

L L

sym L

 − 

 

 − 

 

=  

 − 

 

 

k

KQ = F

(3)

str.2

II. Przykład modelowania układu belkowego

Dana jest belka przedstawiona na rysunku poniżej. Belka o przekroju prostokątnym podparta jest na końcach za pomocą podpór przegubowych, w środku obciążona siłą skupioną P

Model dyskretny składa się z 4 elementów skończonych (5 węzłów), tak więc układ w sumie posiada 10 stopni swobody. Na rysunku zaznaczono również warunki brzegowe dla zadania. Odebranie przemieszczeń liniowych dla węzła 1 oraz 5 (indeksy (1 i 9) oraz siłę w węźle nr 3 (indeks 5)

Na poniższym rysunku zaznaczono poglądowo macierzowy układ równań dla przykładu. Jasnoszare pola w globalnej macierzy sztywności odpowiadają kolejnym lokalnym macierzom sztywności elementu belkowego. Ciemnoszarym kolorem zaznaczono miejsce gdzie wartości ulegają zsumowaniu. Jak wspomniano wcześniej dotyczy to obydwu stopni swobody w węźle. Czarne pola w macierzy sztywności to miejsca w których modyfikacji podlega (zgodnie z zasadami metody kary) element macierzy w którym odebrano przemieszczeniowy stopień swobody. Czarne pole w globalnym wektorze obciążeń odpowiada zadanej sile P.

Uwaga:

Dla utwierdzeń, oprócz przemieszczeń liniowych odbierane są również przemieszczenia kątowe.

W takim przypadku należy to odpowiednio uwzględnić (indeksy parzyste). W wektorze obciążeń parzyste indeksy odpowiadają momentom skupionym w odpowiednich węzłach

(4)

str.3

III. Kod do programu SciLab dla powyższego przykładu

Dane do przykładu:

Moduł Younga E = 200 000MPa Siła P = 200N

Długości poszczególnych elementów: L = 100mm

Parametry przekroju poprzecznego belki b = 5mm, h = 10mm

clear

clc //określenie wartości obciążenia P=200;

//określenie długości elementu (podział jednorodny) L=100;

//określenie własności materiałowych E=2e5;

//obliczenie momentu bezwładności elementów (taki sam dal całej belki) b=5;

h=10;

Iz=b*h^3/12;

//wyznaczenie macierzy sztywności poszczególnych elementów K=E*Iz*[12/L^3, 6/L^2, -12/L^3, 6/L^2;

6/L^2, 4/L, -6/L^2, 2/L;

-12/L^3, -6/L^2, 12/L^3, -6/L^2;

6/L^2, 2/L, -6/L^2, 4/L;];

//Zainicjowanie globalnej macierzy sztywności KG=zeros(10,10);

//Agregacja globalnej macierzy sztywności KG(1:4,1:4)=KG(1:4,1:4)+K;

KG(3:6,3:6)=KG(3:6,3:6)+K;

KG(5:8,5:8)=KG(5:8,5:8)+K;

KG(7:10,7:10)=KG(7:10,7:10)+K;

//Określenie wektora sił węzłowych F=[0; 0; 0; 0; P; 0; 0; 0; 0; 0];

//Określenie przemieszczeniowych warunków brzegowych - metoda kary C=max(KG)*10^4; //Określenie stałej C modelującej sprężynę o dużej sztywności KG(1,1)=KG(1,1)+C;

KG(9,9)=KG(9,9)+C;

//Rozwiązanie układu równań Q=inv(KG)*F;

disp(Q)