AEu = q x (0,l) u(0) = 0 AEu (l)n(l) = S

notatki dowykªadu dla2roku I stopnia kierunku Budowni two

Witold Ce ot

Contents

1 Motywa ja 2

2 Aproksyma ja MES 2

3 Sformuªowanie mo ne zadania prta w stanie jednoosiowym 4

4 Sformuªowania sªabe zadania prta w stanie jednoosiowym 5

5 Metoda Galerkina 7

6 Algorytm MES na przykªadzie zadania 1D 8

7 MES dla sta jonarnego przepªywu iepªa w 2D 16

8 MES dla liniowej spr»ysto± i 24

9 Analiza kratowni za pomo ¡ MES 28

10 MES dla belek 35

11 Zastosowanie MES do zada« dynamiki 38

12 Osza owanie bªdu aproksyma ji MES 43

Naj z± iej obe niestosowan¡ metod¡doobli zaniainteresuj¡ y h projektantawielko± i, taki h jak

ugi ia,siªy przekrojowe,napr»enia, przemiesz zenia jest Metoda Elementów Sko« zony h (MES).

Stanowipodstawalgorytmówu»ywany h wprograma hsªu»¡ y h doprojektowania konstruk jiin-

»ynierski h(np. Robot,Abaqus, Ansys,Catia, AltairHyperWorksiwiele inny h). Poniewa»jednak

zawyniki obli ze«odpowiada projektant korzystaj¡ y z oprogramowania, a nie programista, który

opra owaªsystem, znajomo±¢ podstawMES jestklu zowadlawªa± iwego ibezpie znegostosowania

kodówkomer yjny h.

2 Aproksyma ja MES

W MES stosowana jest aproksyma ja za pomo ¡spe jalny h funk ji bazowy h, zwany h funk jami

ksztaªtu.

Dla funk ji jednej zmiennej

u : [a, b] ∋ x → u(x) ∈ R

^szukamy

u h (x) ≈ u(x)

^{w posta i}

u h (x) = α 1 ϕ 1 (x) + α 2 ϕ 2 (x) + · · · + α n ϕ n (x)

^{, gdzie}

α i

^{to nieznane parametry zwane stopniami swobody, a}

ϕ i (x)

^{s¡funk jamiksztaªtu} deniowanymi za pomo ¡elementów sko« zony h.

1. Dyskretyza ja i liniowe, i¡gªeglobalne funk je ksztaªtu.

1 2 3

1 2 3 4

PSfrag repla ements

ϕ ₁ ϕ ₂ ϕ ₃ ϕ 4

ϕ ^′ ₁ ϕ ^′ ₂ ϕ ^′ ₃ ϕ ^′ ₄

1 1 1 1

1 h 1 h 1 h 1 h

− h ¹

− h ¹

− h ¹

− h ¹

Figure 1: Dyskretyza ja za pomo ¡ 3 elementów oraz wykresy funk ji ksztaªtu Lagrange'a stopnia 1 (s¡

i¡gªe i maj¡ maªe no±niki). Po hodna wsensieklasy znymnie wszdzie istnieje.

2. Deni ja stopni swobody

ϕ i

^{- wielomiany}Lagrange'a,

ϕ i (x k ) = δ ik ⇒ αi = uh(xi)

Je»eli

α _i = u(x _i )

^to

u _h (x)

^jest interpolantem, takim»e

u _h (x _i ) = u(x _i )

^.

3. Pojedyn zy elementsko« zony, algorytmobli ze« elementpoelemen ie

PSfrag repla ements

ϕ 1

ϕ 2

x 1 x 2

Figure 2: Lokalnie ponumerowane funk je ksztaªtu

ϕ 1 (x), ϕ 2 (x)

^elementu sko« zonego o wzªa h

x 1 , x 2

(równie»z zastosowaniem numera jilokalnej).

4. Bª¡daproksyma ji

5. Aproksyma ja Lagrange'a wy»szego stopnia

6. Aproksyma jahierar hi znawy»szego stopnia-

ϕ 1 (x), ϕ 2 (x)

^liniowe,

α 3

minimalizuje

||u − u h ||

7. Funk je ksztaªtu: liniowoniezale»ne, jednozna znierozwi¡zywalne (unisolvent; wektory

ϕ 1 (X), ϕ 2 (X), . . . , ϕ n (X)

^{s¡liniowoniezale»ne}

∀X = (x 1 , x 2 , . . . , x n ) : x i 6= x j

^dla

i 6= j

^).

3 Sformuªowanie mo ne zadania prta w stanie jednoosiowym

WprowadzeniedoMES bdzieprzedstawione naprzykªadziemo»liwieprostego zadaniamodelowego

me haniki jakimjest prt poddany dziaªaniujedynie siªdziaªaj¡ y hw kierunku jego osi.

Równa« opisuj¡ e jednoosiowystan napr»enia w pr ie zyli zadanieroz i¡ganiaalbo± iskania

prta ªatwo wyprowadzi¢ korzystaj¡ z zasady pdu (II zasady dynamiki Newtona)i nastpuj¡ y h

zaªo»e«:

•

^{prt jestw spo zynku (statyka)}

•

przemiesz zeniai odksztaª enia s¡ maªe

•

^{materiaªprta jestspr»ysty (stosujemy prawoHooke'a)}

•

^siªymidzy z¡ste zkowe s¡krótko zasigowe

•

^siªyobjto± iowe (masowe)s¡ i¡gªe

•

^{AE= onst,A - poleprzekroju} poprze znego, E -moduª Younga

•

^{istniej¡ drugiepo hodne} rozwi¡zania wka»dym punk ie od inka

(0, l)

PSfrag repla ements

l[m]

S[N ] q[N/m]

Figure3: Modelowezadanie1D.Prtwjednoosiowymstanienapr»enia. Tylkosiªaosiowamo»eby¢

niezerowa. Pozostaªe siªy przekrojowe (momenty zginaj¡ e, skr aj¡ y i siªy poprze zne) s¡ równe

zero. Przyjmijmy, »e o± "x" ukªadu wspóªrzdny h pokrywa si z osi¡ prta i za zyna na lewym

ko« u.

•

^{prt nieulega} wybo zeniu.

Przytaki hzaªo»enia h,dlaprtajaknarys. 3przemiesz zenia

u(x)

^{przekrojówprtawkierunku}

jegoosimo»naobli zy¢rozwi¡zuj¡ nastpuj¡ ezagadnieniebrzegowe(równanieró»ni zkoweiwarunki

brzegowe)







−AEu ^′′ = q ∀x ∈ (0, l) u(0) = 0

AEu ^′ (l)n(l) = S

(1)

gdzie

n(l) = 1, n(0) = −1

^{, ob i¡»enia}

q, S

^{maj¡ znakizgodne z przyjtym ukªadem} wspóªrzdny h.

Jest to sformuªowanie mo ne, którego e h¡ jest speªnianie w ka»dym punk ie obszaru (dziedziny

funk ji

u

^{zyli od inka}

[0, l]

^{) pewnego warunku.}

4 Sformuªowania sªabe zadania prta w stanie jednoosiowym

4.1 Sformuªowanie waria yjne

Metoda elementów sko« zony h wymaga t.zw. sformuªowania sªabego. Poni»ej przedstawiona jest

kon ep ja przeksztaª enia sformuªowania mo nego w sªabe z pomini iem sz zegóªów matematy-

zny h.

Pomnó»my równanieró»ni zkowe (1)

i

obustronnie przezdowoln¡ i¡gª¡funk j

v(x)

− vAEu ^′′ = vq ∀v ∈ V

⁽²⁾

Otrzymali±my w ten sposób nowe równanie albo ina zej pewn¡ funk j zmiennej

x

przedstawion¡

nadwa sposoby: zapomo ¡drugiej po hodnejniewiadomej funk ji

u(x)

⁽

−vAEu ^′′

^{)orazzapomo ¡}

ob i¡»enia (

vq

^{). Zatem aªki ozna zone w grani a h od}

0

^do

l

^{z wyra»e« po lewej i prawej stronie}

równania(2)bd¡ równe zyli

− Z l

0

vAEu ^′′ dx = Z l

0

vq dx ∀v ∈ V

⁽³⁾

Po zastosowaniu aªkowania przez z± i do aªkipolewej stronie otrzymamy

Z l 0

v ^′ AEu ^′ dx − [v(l)AEu ^′ (l) − v(0)AEu ^′ (0)] = Z l

0

vq dx ∀v ∈ V

⁽⁴⁾

Funk je

v

^{zwane funk jami testowymi byªy do tej pory dowolnymi funk jami i¡gªymi (elementy}

V

^{). Je»eli} ograni zymy si do funk ji które maj¡ warto±¢ 0 dla

x = 0

przeksztaª enie ni nie stra i na ogólno± i, gdy» ty h spe jalny h funk ji testowy h jest w i¡» niesko« zenie wiele. Natomiast

warunekzerowaniasifunk jitestowejwmiejs uwktórymrozwi¡zaniemaznan¡warto±¢nabrzegu

(speªnia warunek kinematy zny na podporze o wzapisie matematy znym zazna zamy pisz¡ , »e

v

nale»ydoprzestrzeni

V 0

^{)uprosz zawistotnysposób}sformuªowanie,którepouwzgldnieniuwarunku staty znego (

AEu ^′ (l) = S

^{) przyjmuje posta¢}

Znale¹¢

u(x), u(0) = 0

^{takie, »e}

Z l

0

v ^′ AEu ^′ dx = Z l

0

vq dx + v(l)S ∀v ∈ V 0

⁽⁵⁾

Zauwa»my,»ekinematy znywarunekbrzegowyzostaªprzepisany(nazywasiina zejpodstawowym),

a staty zny jest uwzgldniony po prawej stronie (st¡d nazwa naturalny). Sformuªowanie (5) jest

dobrze znan¡ z me haniki zasad¡ pra wirtualny h. Funk je

v

^{ozna zaj¡} przemiesz zenia wirtu- alne ( i¡gªe i zeruj¡ e si na podpora h), lewa strona tej zasady to pra a siª wewntrzny h na

przemiesz zenia h wirtualny h a wyra»enie po prawej stronie jest pra ¡ siª zewntrzny h na ty h

przemiesz zenia hwirtualny h. Jest to sformuªowanie sªabe, gdy» zasada(5) niezakªadaspeªniania

pewny hwarunkówwpunkta h. Jegopodstawowezaletytoni»szystopie«po hodnejzposzukiwanej

funk ji i dobre podstawy matematy zne. Ponadto, okazuje si, »e zasada obowi¡zuje dla zmiennej

sztywno± i

AE

^oraz nie i¡gªego ob i¡»enia

q

^.

4.2 Problem minimaliza ji

Zasad pra wirtualny h mo»na przeksztaª i¢ do równowa»nej zasady minimum aªkowitej energii

poten jalnej (

J : V 0 ∋ v → J(v) ∈ R

^{)1. Dla prta poddanego} przemiesz zeniom wirtualnym

v

wynikaj¡ a ztego aªkowita energiapoten jalna

J(v)

^{wyra»a si wzorem}

J(v) = 1

2 Z l

0

v ^′ AEv ^′ dx − Z l

0

vq dx − v(l)S

⁽⁶⁾

Przemiesz zenia osiprta

u(x)

^s¡tymiprzemiesz zeniami wirtualnymi(

v(x), v(0) = 0

^{) dlaktóry h}

funk jonaª

J

^{osi¡gaminimum zyli}

J(u) = min _V 0 J(v)

⁽⁷⁾

Podsumowuj¡ , przemiesz zenia osi prta mo»na obli zy¢ rozwi¡zuj¡ zadanie (1) lub (5) lub (7).

Wszystkie maj¡ jednozna zne rozwi¡zanie przy odpowiedni hwarunka hbrzegowy h.

5 Metoda Galerkina

Zewzgldunanajwiksz¡ogólno±¢zastosujmy zasadpra wirtualny h(5) zylisformuªowaniewari-

a yjne. Przybli»one rozwi¡zanie (

u h (x)

^{) takiego zadaniamo»na otrzyma¢ metod¡ Galerkina, który}

ponad 100 lat temu zaproponowaª aby rozwi¡zanie przybli»one poszukiwa¢ w posta i kombina ji

1

Formalnie dzikitemu, »e funk jonaªna pra  siª zewntrzny h na przemiesz zenia h wirtualny h (lewastrona

(5)jestsymetry znyidodatniookre±lony.

liniowej pewny h znany h funk ji bazowy h ( i¡gªy h i liniowoniezale»ny h)

ϕ 1 (x), ϕ 2 (x), ..., ϕ n (x)

zyli

u(x) ≈ u h (x) = α 1 ϕ 1 (x) + α 2 ϕ 2 (x) + · · · + α n ϕ n (x)

⁽⁸⁾

gdzie

α i

^{s¡ nieznanymi} wspóª zynnikami (stopniami swobody). Funk je bazowe

ϕ i

^{s¡ równie»}

stosowane jako przemiesz zenia wirtualne 2

. Przykªadowo przyjmuj¡

v = ϕ 1

^{mamy pierwsze rów-}

naniealgebrai zne

Z l 0

ϕ ^′ ₁ AE(α ₁ ϕ ^′ ₁ + α ₂ ϕ ^′ ₂ + · · · + α _n ϕ ^′ _n ) dx = Z l

0

ϕ ₁ q dx + ϕ ₁ (l)S

⁽⁹⁾

Zapisuj¡ równanie (5) dlakolejny h funk ji bazowy h traktowany h jako funk je testowe otrzymu-

jemy nastpuj¡ y symetry zny ukªad

n

^równa« algebrai zny hz

n

niewiadomymi

Ga = f

⁽¹⁰⁾

gdzie

a

^{jestkolumnowym wektorem} niewiadomy hstopni swobody

α 1 , α 2 , . . . , α n

^,

G ij = Z l

0

ϕ ^′ _i AEϕ ^′ _j dx f i = Z l

0

ϕ i q dx + ϕ i (l)S i, j = 1, 2, ..., n

⁽¹¹⁾

Rozwi¡zanie powy»szego ukªadu równa« algebrai zny h pozwala obli zy¢ stopnieswobody, a zetem

równie»funk j

u h (x)

^{napodstawie wzoru (8).}

Takasamaaproksyma jajakwmetodzieGalerkinadlazasadyminimumfunk jonaªu(7)prowadzi

do metody Ritza, która daje ten sam ukªad równa« algebrai zny h (10) o metoda Galerkina.

Poprawne sformuªowanie zadania i i¡gªe funk je bazowe daj¡ gwaran j zbie»no± i metody zyli

zmniejszanie bªduwraz ze zwikszaniem li zby stopni swobody.

6 Algorytm MES na przykªadzie zadania 1D

Metoda elementów sko« zony h pozwala uzyskiwa¢ przybli»one, ale z kontrolowan¡ dokªadno± i¡,

rozwi¡zania zagadnie« brzegowy h modeluj¡ y h za howanie si rze zywisty h konstruk ji. Jest

metod¡Galerkina,wktórejjakofunk jibazowy hu»ywasipoznaneju»w ze±niejfunk jeksztaªtu.

W tym rozdzialeprzedstawiony jestalgorytm MES naprzykªadzie zadaniaprta.

1. Sformuªowanie zadania

Skorzystamy ze sformuªowania waria yjnego - zasady pra wirtualny h (5). Wiemy, na pod-

stawieanalizyfunk jonalnej,»ejesttodobrzepostawionezadanie,maj¡ ejednozna znerozwi¡-

zanie. Przyjmijmy,»e: L=2 m,AE=1000kN, q=10kN/m,S=-3kN. Dlawszystki htypowy h

zada«sformuªowania s¡ dobrzeznane i u»ytkownik programówMES niemusii h podawa¢.

2. Dyskretyza ja - genera ja siatki

Podzielmy od inek

[0, l]

^{, w którym} poszukiwane jest rozwi¡zanie na mniejsze od inki, bd¡ e w tym przypadku elementami sko« zonymi. Jest to t.zw. genera ja siatki MES. Zaªó»my, »e

2

Istnieje wersjametody,wktórejfunk jamitestowymi(przemiesz zeniamiwirtualnymi)s¡innefunk jeni»przyj-

mowanedoaproksyma jirozwi¡zania.

zastosujemy3elementysko« zone (rys. 4). Generatorsiatkidostar zainforma jiowspóªrzd-

ny h wzªów (

x 1 = 0, x 2 = ¹ ₃ l, x 3 = ² ₃ l, x 4 = l

^{) i o tym jakie wzªy deniuj¡} posz zególne elementy sko« zone. Istotna jest kolejno±¢ wzªów przyporz¡dkowany h elementom. Zaªó»my

nastpuj¡ e rela je: element1 ma wzªy (1,2),element 2 mawzªy (2,3) i element 3ma wzªy

(3,4). Wprowadzony podziaª zyli dyskretyza ja umo»liwia zdeniowanie funk ji ksztaªtu. S¡

1 2 3

1 2 3 4

PSfrag repla ements

ϕ ₁ ϕ 2

ϕ ₃ ϕ ₄

ϕ ^′ ₁ ϕ ^′ ₂ ϕ ^′ ₃ ϕ ^′ ₄

1 1 1 1

1 h 1 h 1 h 1 h

− h ¹

− h ¹

− h ¹

− h ¹

Figure4: Modelowezadanie1D.Dyskretyza jazapomo ¡3elementówsko« zony h zwprowadzon¡

numera j¡ elementów i wzªów oraz wykresy funk ji ksztaªtu i i h po hodny h. Element 1 ma

wzªy (1,2), element 2 ma wzªy (2,3), a element 3 ma wzªy (3,4). Wspóªrzdne wzªów s¡ znane.

Ozna zmy je odpowiednio

x 1 , x 2 , x 3 , x 4

^. Przyjmijmy, »e dªugo± i elementóws¡ jednakowe i wynosz¡

h = ₃ ^l

^{. Wtedy}

x ₁ = 0, x ₂ = ¹ ₃ l, x ₃ = ² ₃ l, x ₄ = l

^.

tofunk je i¡gªe(dzikitemumetodanumery znajeststabilna zylimamygwaran j,»ezwik-

szaj¡ li zbelementówsko« zony hzmniejszamybª¡daproksyma jiMES),przyporz¡dkowane

wzªomwtakisposób,»eka»dazfunk jiksztaªtuprzyjmujewarto±¢1wwybranymw¹leoraz0

wpozostaªy hwzªa h. Przyjmijmy,»etakjaknarys. 4bd¡tofunk jeod inkowoliniowe(i h

wykresami s¡ ªamane) 4

. Funk je ksztaªtu bd¡ funk jami bazowymi aproksyma ji rozwi¡za-

nia,funk jamitestowymiwmetodzieGalerkinaorazfunk jamiaproksymuj¡ ymibrzegobszaru

(w 1D zawsze, w 2D i 3D w wielu przypadka h aproksyma ja brzegu jest dokªadnym opisem

geometriiobszaru).

Aproksyma ja MES wyra»a si wzorem (8). Dziki opisanym powy»ej wªa± iwo± iom funk ji

ksztaªtuªatwozauwa»y¢, »e

u h (x 1 ) = α 1 ϕ 1 (x 1 ) + α 2 ϕ 2 (x 1 ) + · · · + α n ϕ n (x 1 ) = α 1 · 1 + α 2 · 0 + · · · + α n · 0 = α 1

⁽¹²⁾

orazogólnie

u h (x i ) = α i , i = 1, 2, . . . , n

⁽¹³⁾

3

W prakty epotrzebny hjestzna zniewi ej elementówdouzyskaniaodpowiedniejdokªadno± i.

4

Klasy znepo hodnefunk jiksztaªtunieistniej¡wwzªa h(prawostronnanierównasilewostronnej)alezarówno

po hodnejaki aªkiwzasadziewaria yjnej(5)s¡rozumianewogólniejszymsensieidlategozmatematy znegopunktu

widzeniawszystkojestpoprawnie. Zpunktuwidzeniame hanikiprzemiesz zeniamog¡by¢niegªadkie(mie¢skokow¡

zmianpo hodnej) wmiejs uzmianyparametrówmateriaªowy h o mamiejs enp. wkompozyta h.

Wprowadzaj¡ ozna zenie

u h (x i ) = u i

aproksyma j MES przyjmiemy w posta i kombina ji funk jiksztaªtuzestopniamiswobodybd¡ ymiwarto± iami(narazienieznanymi)rozwi¡zania

przybli»onegow wzªa h

u h (x) = u 1 ϕ 1 (x) + u 2 ϕ 2 (x) + · · · + u n ϕ n (x)

⁽¹⁴⁾

Podstawow¡ zalet¡ funk ji ksztaªtu jako funk ji bazowy h jestto »e przyjmuj¡ warto± i nieze-

rowe tylko na niewielki h podobszara h (w elementa h zwi¡zany h z wybranym wzªem). W

prakty e, przy du»ej li zbie elementów sko« zony h, podobszary te s¡ bardzo maªe w porów-

naniuz aªym obszarem.

W rozwa»anym przykªadzie

n = 4

^{i mówimy, »e globalna li zba stopni swobody (LSS) jest}

równa 4.

3. Agrega ja

Poniewa»metodaelementówsko« zony h jestmetod¡Galerkinaukªad równa«algebrai zny h

pozwalaj¡ yobli zy¢stopnieswobodymaposta¢tak¡jak wrozdziale5alez dodatkow¡ enn¡

zalet¡jak¡jest pasmowo±¢ (du»o wspóª zynników zerowy h).

Przykªadowoje»elijakofunk jitestoweju»yjemypierwszejfunk jiksztaªtuotrzymamypierwsze

równaniealgebrai zne

Z l 0

ϕ ^′ ₁ AE(u 1 ϕ 1 + u 2 ϕ 2 + u 3 ϕ 3 + u 4 ϕ 4 ) ^′ dx = Z l

0

ϕ 1 qdx + ϕ 1 (l)P

⁽¹⁵⁾

Poniewa»

ϕ ^′ ₁ ϕ ^′ ₃ = ϕ ^′ ₁ ϕ ^′ ₄ = 0

^,a

ϕ ^′ ₁ ϕ ^′ ₁

ⁱ

ϕ ^′ ₁ ϕ ^′ ₄

^{s¡niezerowetylkowelemen iepierwszym,powy»sze}

równaniemaposta¢

u ₁ Z

e 1

ϕ ^′ ₁ AEϕ ₁ ^′ dx + u ₂ Z

e 1

ϕ ^′ ₁ AEϕ ₂ ^′ dx + u ₃ · 0 + u ₄ · 0 = Z

e 1

ϕ ₁ qdx + 0 · S

⁽¹⁶⁾

gdzie

e 1

^{ozna za pierwszyelement(od inek}

[0, h]

^{). Nastpnerównia}otrzymujemyanalogi znie imo»na je zapisa¢ma ierzowo

AE























 R

e1 ϕ ^′ ₁ ϕ 1 ′ R

e1 ϕ ^′ ₁ ϕ 2 ′ 0 0

R

e1 ϕ ^′ ₂ ϕ 1 ′ R

e1 ϕ ^′ ₂ ϕ 2 ′ + R

e2 ϕ ^′ ₂ ϕ 2 ′ R

e2 ϕ ^′ ₂ ϕ 3 ′ 0

0 R

e2 ϕ ^′ ₃ ϕ 2 ′ R

e2 ϕ ^′ ₃ ϕ 3 ′ + R

e3 ϕ ^′ ₃ ϕ 4 ′ R

e3 ϕ ^′ ₃ ϕ 4 ′

0 0 R

e3 ϕ ^′ 4 ϕ 3 ′ R

e3 ϕ ^′ 4 ϕ 4 ′











































 u 1

u 2

u 3

u 4





















=

















 R

e1 ϕ 1 q R

e1 ϕ 2 q + R

e2 ϕ 2 q R

e2 ϕ 3 q + R

e3 ϕ 3 q R

e3 ϕ 4 q

















 +















 0

0

0

S

















(17)

gdzie dla skró enia zapisu pominito pod aªkami symbol

dx

^{oraz przyjto, »e}

AE =

^{onst, a}

kolorami ozna zono wkªad kolejny h elementów sko« zony h do (globalnego) ukªadu równa«

algebrai zny h. Podma ierze2x2ozna zoneposz zególnymikoloraminazywaj¡ sima ierzami

(sztywno± i)elementów. Analogi zniepodma ierzekolumnowe2x1w hodz¡ ewskªadwektora

prawy h stron nazywaj¡ si wektorami(ob i¡»enia)elementów. Dla ka»degoelementu korzys-

tamy wzwi¡zkuztymjedyniezdwó hfunk jiksztaªtu,gdy» tylkote, któres¡zwi¡zanezjego

wzªami(rys. 5)niezeruj¡siwelemen ie. Dla ty hfunk jistosujemynumera jlokaln¡: 1,2.

ZatemLSS elementuw tym zadaniui przy tej dyskretyza ji jest równa 2.

PSfrag repla ements

ϕ ₁ ϕ ₂

Figure 5: Modelowe zadanie 1D. Lokalnie ponumerowane funk je ksztaªtu

ϕ 1 (x), ϕ 2 (x)

^elementu

sko« zonego o wzªa h

x 1 , x 2

^{(równie» z}zastosowaniem numera ji lokalnej).

Stosuj¡ lokaln¡numera j(niebdzieonadalejspe jalnieozna zana),wspóª zynnikima ierzy

(

K ^e

^{)i wektora (}

P ^e

^)elementu

e

^{obli zasi ze wzorów}

K _ij ^e =

Z

e

ϕ ^′ _i AEϕ ^′ _j dx, P _i ^e = Z

e

ϕ i q dx, i, j = 1, 2

⁽¹⁸⁾

Wzory te mo»na uzyska¢ bezpo±rednio, podstawiaj¡ funk je ksztaªtu za funk je testowe i

bazowe do aªek po lewej i prawej stronie sformuªowania waria yjnego. Skªadnik bez aªki

(

v(l)S

^{) nie jest zwi¡zany z »adnym elementem i jest} uwzgldniany ina zej ni» skªadniki z aªkami.

Šatwoobli zy¢ma ierzsztywno± i dlaelementuodªugo± i

h

(niezale»nieodjegopoªo»eniana osi

x

^{) oraz wektor ob i¡»enia dlastaªegoob i¡»enia}

q

K ^e = AE h

 1 −1

−1 1



, f ^e = qh 2

 1 1



(19)

Ma ierzwspóª zynników (

K

^{)iwektor prawy h stron(}

f

⁾algebrai znego ukªadurówna« nazy- waj¡ siglobaln¡ma ierz¡(sztywno± i)iglobalnymwektorem(ob i¡»enia). I hwspóª zynniki

najwygodniej obli za si metod¡ agrega ji, element po elemen ie 5

. Agrega j w nastpuj¡ y

sposób.

•

Przyjmujemy zerow¡globaln¡ma ierzsztywno± i

K _nxn

^{izerowyglobalnywektorob i¡»e-}

nia

f _nx1

^{. Poniewa»}

n = 4

K =









0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0









, f =







 0 0 0 0









(20)

•

^{Kolejno, dla ka»dego elementu obli zamy ma ierz i wektor elementowy6, a nastpnie do-}

dajemy i h wspóª zynniki do ma ierzy i wektora globalnego w miejs a wynikaj¡ e z nu-

merówglobalny h stopni swobody elementu. Numery te s¡ adresami(okre±laj¡ kolumny

i wiersze) w ma ierzy i wektorze globalnym pod które nale»y "dostar zy¢" posz zególne

wspóª zynniki ma ierzy iwektora elementu.

5

Przykªad 1Djestsz zególny,gdy» sformuªowaniewaria yjnezawieraskªadnikbez aªki, któryuwzgldniasina

poziomieglobalnym, zylibezagrega ji.

6

W rozwa»anymprzykªadziema ierzeiwektorywszystki helementóws¡takiesame(patrzwzory19).

K = AE h









1 −1 0 0

−1 1 0 0 0 0 0 0 0 0 0 0









, f = qh 2







 1 1 0 0









(21)

(b) Element 2 mawzªy (2,3). Po agrega jijego ma ierzyi wektora

K = AE h









1 −1 0 0

−1 2 −1 0 0 −1 1 0

0 0 0 0









, f = qh 2







 1 2 1 0









(22)

( ) Element 3 mawzªy (3,4). Po agrega jijego ma ierzyi wektora

K = AE h









1 −1 0 0

−1 2 −1 0

0 −1 2 −1

0 0 −1 1









, f = qh 2







 1 2 2 1









(23)

•

^{W nietypowym zadaniu 1D nale»y jesz ze doda¢ wektor z siª¡ skupion¡}

S

^{przyªo»on¡ w}

4-tym w¹le (dlategowystpuje w4-tym wierszu), a zatem

K = AE h









1 −1 0 0

−1 2 −1 0

0 −1 2 −1

0 0 −1 1









, f = qh 2







 1 2 2 1







 +







 0 0 0 S









(24)

Otrzymali±myw ten sposób ukªad równa« algebrai zny h

AE h









1 −1 0 0

−1 2 −1 0

0 −1 2 −1

0 0 −1 1















 u 1

u 2

u ₃ u 4









= qh 2







 1 2 2 1







 +







 0 0 0 S









(25)

Przyjmuj¡ ozna zenia







 u 1

u 2

u 3

u 4









= u, qh 2







 1 2 2 1









= P ,







 0 0 0 S









= W

⁽²⁶⁾

zyli wektory, odpowiednio: stopni swobody, ob i¡»enia i siª wzªowy h mo»na zapisa¢ ukªad

równa«algebrai zny h w zwarty sposób

Ku = P + W

⁽²⁷⁾

0 0 0 0 0

0 0 0 0 0

0 0 0 0 0

0 0 0 0 -3

1500 -1500 0 0 3.3333

-1500 1500 0 0 3.3333

0 0 0 0 0

0 0 0 0 -3

1500 -1500 0 0 3.3333

-1500 3000 -1500 0 6.6667

0 -1500 1500 0 3.3333

0 0 0 0 -3

1500 -1500 0 0 3.3333

-1500 3000 -1500 0 6.6667

0 -1500 3000 -1500 6.6667

0 0 -1500 1500 0.33333

Figure6: Modelowe zadanie1D. Przebieg agrega ji ma ierzy sztywno± i (4x4)i wektora ob i¡»enia

(ostatniakolumna). Wektor siªwzªowy h wstawionona po z¡tku agrega ji.

Po przyj iu parametrów zdeniowany h na po z¡tku rozdziaªu 6, ma ierz i wektor elementu

s¡równe

K ^e =

 1500 −1500

−1500 1500



, P ^e =  3.333 3.333



(28)

iagrega ja przebiega wsposób pokazany narys. 6.

4. Uwzgldnienie warunków kinematy zny h i rozwi¡zanie ukªadu równa«

Ma ierz ukªadu równa« algebrai zny h (27) jest osobliwa. Wynika to z nie uwzgldnienia do

tej pory kinematy znego warunku brzegowego. Z me haniki wiemy, »e bez taki h warunków

(podpór) ukªad jest geometry znie zmienny i analiza staty zna nie umo»liwia obli zy¢ jed-

nozna zny h przemiesz ze«. Zmatematy znegopunktu widzenia warunekkinematy zny, zyli

u(0) = 0

^{implikuje,»epierwszafunk jaksztaªtuniejestwªa± iw¡funk j¡testow¡}(przemiesz ze- niem wirtualnym), a zatem pierwsze równanie w uzyskanym metod¡ agrega ji ukªadzie (27)

jest"nielegalne". Zostaªoono jednak utworzone elowo, dlauprosz zenia algorytmuagrega ji

ale teraz nale»y go wyeliminowa¢. Poniewa» znamy przemiesz zenie w pierwszym w¹le zyli

znamy pierwszy stopie« swobody

u 1

^{mo»na pomno»y¢ pierwsz¡ kolumn ma ierzy sztywno± i}

przez

u 1

^{(gdyby byª znany stopie« swobody}

u k

^{to kolumn} "k") i przenie±¢ otrzymany wek- tor na praw¡ stron

7

. W tym przykªadzie

u ₁ = 0

^zatem przenoszony wektor jest zerowy, ale niezawsze tak musiby¢. Nastpnierezygnujemy z pierwszego (w ogólno± ik-tego) równania.

Otrzymujemy w ten sposób ukªad 3równa« z 3niewiadomymi

AE h





2 −1 0

−1 2 −1 0 −1 1







 u 2

u 3

u 4



 = qh 2



 2 2 1



 +



 0 0 S





⁽²⁹⁾

Numery znie, uwzgldnienie kinematy zny h warunków brzegowego wygodniej jest wykona¢

bez zmiany wymiarów ma ierzy i wektora ukªadu równa«. Po przeniesieniu, jak poprzednio,

znany h wyra»e« na prawe strony nale»y wyzerowa¢ k-t¡ kolumn (przeniesiona na praw¡

stron) i k-ty wiersz ("nielegalne" równanie) ma ierzy sztywno± i a nastpnie wstawi¢ li zb

1

^{na przek¡tn¡ oraz znany stopie« swobody na praw¡ stron k-tego równania (warunek}

u k =

7

Jest toprzenoszeniewka»dymrównaniuznanegowra»enianapraw¡stron.

u(x k )

^{). W rozwa»anym} przykªadzie otrzymamy

AE h









1 0 0 0

0 2 −1 0

0 −1 2 −1

0 0 −1 1















 u 1

u 2

u 3

u 4









= qh 2







 0 2 2 1







 +







 0 0 0 S









(30)

Ukªad równa« algebrai zny h, w posta i (29) albo (30) jest symetry zny, pasmowy, dodatnio

okre±lony i ma jednozna zne rozwi¡zanie, które mo»na obli zy¢ jedn¡ z bezpo±redni h albo

itera yjny hmetodnumery zny h. Dlaprzyjty hkonkretny h dany hli zbowy hotrzymamy

u 1 = 0, u 2 = 0.009111, u 3 = 0.01378, u 4 = 0.01400 [m]

⁽³¹⁾

5. Wyniki obli ze«

Znaj¡ stopnie swobody mo»emy wró i¢ do aproksyma ji za pomo ¡ funk ji ksztaªtu. Na-

jlepiejstosowa¢ j¡element poelemen ietak jaktobyªoprzy agrega jiukªadu równa«. Zatem,

ka»dy element sko« zony rozwa»amy indywidualnie, stosuj¡ numera j lokaln¡ (rys. 5).

Wiedz¡ , na podstawie rela jielement-wzªy jakie s¡ numery wzªów pobieramy wspóªrzdne

ty h wzªów (w lokalnej numera ji ozna zamy je

x ₁ , x ₂

^), od zytujemy z globalnego wektora stopni swobody stopnie swobody elementu (ozna zaj¡ jako

u 1 , u 2

^{) i zapisujemy} rozwi¡zanie MESw rozwa»anym elemen ie sko« zonym

u ^e _h (x) = u 1 ϕ 1 (x) + u 2 ϕ 2 (x) = u 1

x − x 2

x ₁ − x ₂ + u 2

x − x 1

x ₂ − x ₁

⁽³²⁾

Przykªadowo, dla elementu 3, stosuj¡ numera j lokaln¡

x 1 = 1.333, x 2 = 2.000, u 1 = 0.01378, u 2 = 0.01400

^{, a zatem}

u ⁽³⁾ _h = x/3000 + 1/75

⁽³³⁾

Mo»emy teraz obli zy¢ siªosiow¡ wka»dym elemen ie korzystaj¡ ze wzoru

N = AEu ^′

^{, zyli}

dla rozwi¡zania przybli»onego (32), z lokanie ponumerowanymi stopniami swobody

u ₁ , u ₂

^,

za hodzi

N _h ^e (x) = AE[u 1 ϕ ^′ 1 (x) + u 2 ϕ ^′ 2 (x)] = AE u 2 − u 1

h

⁽³⁴⁾

Wsz zególno± i, dlaelementu 3:

N _h ⁽³⁾ = AE 0.01400 − 0.01378

2/3 = 0.3333 kN

⁽³⁵⁾

Takuzyskane rozwi¡zaniei siªaosiowa,narysowane elementpoelemen ie,s¡przedstawionena

rys. 7.

6. Bª¡d rozwi¡zania

Po uzyskaniu rozwi¡za« wa»n¡ rze z¡jest znajomo±¢ jego dokªadno± i. Bªdy przemiesz ze« i

siªyosiowej deniuje si nastpuj¡ o

e u (x) = u(x) − u h (x), e N (x) = N(x) − N h (x)

⁽³⁶⁾

0 0.5 1 1.5 2 x

-2 0 2 4 6 8 10 12 14 16

u h 10 ^-3

0 0.5 1 1.5 2

x -4

-2 0 2 4 6 8 10 12 14 16 18

N h

Figure7: Modelowezadanie1D.Wykresrozwi¡zaniaMESiobli zonejnajegopodstawiesiªyosiowej.

0 0.5 1 1.5 2

x -2

0 2 4 6 8 10 12 14 16

u h

10 ^-3

0 0.5 1 1.5 2

x -4

-2 0 2 4 6 8 10 12 14 16 18

N h

Figure8: Modelowezadanie1D.Wykresrozwi¡zaniaMESiobli zonejnajegopodstawiesiªyosiowej

(liniegrube), dokªadne rozwi¡zanie isiªa osiowa(linie kreskowane), bªdy rozwi¡zaniai siªy osiowej

(linie ienkie).

Przykªadowe zadanie, które rozwi¡zujemy zapomo ¡MES jest bardzo prostei ªatwo obli zy¢

dokªadne

u(x), N(x)

^{. I h wykresy pokazano lini¡ przerywan¡ na rys 8.} Zamiesz zono tam równie»wykresybªdów

e u (x), e N (x)

^{. W}rze zywisto± iMESstosujemywtedy,kiedynieznamy rozwi¡zaniadokªadnego. Dlatego wprakty emo»liwejestjedynie osza owanie bªdu obli ze«.

Metody osza owania bªdubez znajomo± i rozwi¡zaniadokªadnego bd¡ przedstawione wdal-

szej z± i wykªadu.

Wartozauwa»y¢, »e przemiesz zenias¡ i¡gªeidokªadnewwzªa hsiatkiMES, natomiastsiªa

osiowa jest od inkowo i¡gªai dokªadna w ±rodka h elementów. Dokªadne przemiesz zenia w

wzªa h wynikaj¡ ze sz zególny h e h zada« 1D. W ogólno± i obserwujemy, »e rozwi¡zanie

MESmanajmniejszybª¡d wwzªa ha napr»eniawpewny h punkta hwewn¡trz elementów.

Bª¡d napr»e« (siª) jest na ogóª wikszy od bªdu przemiesz ze«. Czsto, a zwªasz za dla

ukªadówprtowy hstosujesispe jalneopra owaniewyników,którepozwalaobli zy¢dokªadne

warto± i siª przekrojowy h ukªadów prtowy h w wzªa h siatki MES (zamiast skokowej i h

W przedstawionym powy»ej algorytmie mo»na wyró»ni¢ trzy harakterysty zne etapy obli ze«:

prepro esor (punkty 1,2) -sformuªowaniezadania, przyj ie dany h,genera ja siatkiMES istruk-

tury dany h, pro esor (punkty 3,4) - utworzenie i rozwi¡zanie ukªadu równa« algebrai zny h, za-

pamitaniewyniku, postpro esor(punkty5,6)-przedstawieniewyników(przemiesz ze«)obli zenie

po hodny h (napr»e«),osza owanie bªdu obli ze«.

Algorytm MES zostaª zilustrowany prostym zadaniem me haniki. Jednak interpreta ja zy zna

nie miaªa »adnego zna zenia. Istotny byª typ rozwa»anego równania ró»ni zkowego (

−AEu ^′′ = q

^),

któremo»nazast¡pi¢np.,bardzopodobnym

−AkT ^′′ = f

^,modeluj¡ ymsta jonarnyprzepªyw iepªa (

T (x)

^{ozna zarozkªad}temperatury,

k

^jestwspóª zynnikiemprzewodzenia iepª¡,a

f

intensywno± i¡

¹ródªa iepªa)ialgorytmmetodypozostanieidenty zny. Dlatego,MESjestmetod¡aproksyma ji

rozwi¡za« równa« ró»ni zkowy h, a równania te s¡ matematy znymi modelami zjawisk

zy zny h. Warto o tympamita¢, zwªasz za przy o enie wiarygodno± iwynikówobli ze«.

7 MES dla sta jonarnego przepªywu iepªa w 2D

W rozdziale omówione jest zastosowanie metody elementów sko« zony h do rozwi¡zania równania

Lapla e'a w obszarze 2D. Równanie to mo»e opisywa¢ przepªyw iepªa w iaªa h staªy h, gdzie

zjawiskotoodbywasinazasadzieprzewodzeniaina zejdyfuzji 8

. Dodatkowozaªo»ymy,»e mamydo

zynieniaz pro esem ustalonym, a wi sytua j¡, gdy temperatura niezmienia si w zasie. Zale»y

jedynie od poªo»enia (wspóªrzdny h przestrzenny h) punktu i dla uprosz zenia tylko od dwó h

wspóªrzdny h. Bdziemy zatem poszukiwa¢funk ji dwó h zmienny h o warto± ia h skalarny h.

Równaniami opisuj¡ ymi rozkªad temperatury

T (x, y)

^przy sta jonarnym przepªywie iepªa w iele staªym s¡: zasada za howania energii dlaukªadów termodynami zny h

div

q = Q

⁽³⁷⁾

oraz prawo Fouriera jako najprostsze prawo konstytutywne, które dla materiaªu izotropowego ma

posta¢

q = −k∇T

⁽³⁸⁾

gdzie div

q = ∂q x

∂x + ∂q y

∂y

^,

∇T =  ∂T

∂x , ∂T

∂y



,

q

^jest strumieniem iepªa

[ _m ^W 2 ]

^{, ± i±lej nat»eniem}

strumienia, którego warto±¢

q = lim A,t→0

Q

At

^{(ilo±¢ iepªa} przepªywaj¡ ego przez dan¡ powierz hni

podzielonaprzezpoletejpowierz hni i zas przepªywu),

Q

^{ozna za} intensywno±¢ ¹ródªa iepªa

[ _m ^W 3 ]

^,

któramo»e by¢ dodatnia (ogrzewanie) alboujemna ( hªodzenie),

k

^{jest parametrem} materiaªowym zwanym wspóª zynnikiem przewodzenia iepªa (przykªadowo k=0.8

W

m deg

^{dla betonu, k=50}

W m deg

dla stali). W ogólno± i, np. w iaªa h anizotropowy h, wektor

q

^{nie jest} wspóªliniowy z gradien- tem temperatury i wªa± iwo± i materiaªu s¡ opisywane za pomo ¡ ma ierzy, która w najprostszym

przypadku,dla przyjtegotu materiaªujednorodnego izotropowego (

q || ∇T

^{)ma posta¢}

K =  k 0 0 k



(39)

8

Inneformyprzepªywu iepªa,któres¡opisywaneinnymirównaniamitokonwek jaipromieniowanie.

muªowanie mo ne zadania sta jonarnego przepªywu iepªa w iaªa h staªy h: znale¹¢ temperatur

T(x,y) w obszarze

Ω

^{(rys. 9), tak¡ »e}



 



 



−k∆T = Q

^w

Ω T = ˆ T

^na

∂Ω T

k ^∂T _∂n = q n

^na

∂Ω q

(40)

gdzie

∆T = ∂ ² T

∂x ² + ∂ ² T

∂y ² , ∂T

∂n = ∂T _x

∂x n _x + ∂T _y

∂y n _y , n _x , n _y

^{s¡ skªadowymi wersora (}

pn ² _x + n ² _y = 1

⁾

prostopadªegodobrzegu,skierowanegonazewn¡trz,

q n

^{jestwarto± i¡skªadowejnormalnejstrumienia}

iepªana brzegu, maj¡ ¡znak dodatni je»eli iepªo dopªywa do iaªa (analogi zniejak dla

Q

^{). Aby}

PSfrag repla ements

T ˆ

Q ˆ

q n

ˆ q _n

Figure 9: Ozna zenia stosowane w zadaniu sta jonarnego przewodzenia iepªa.

Q

^{- ¹ródªo iepªa}

wewn¡trz iaªa (ma warto±¢ dodatni¡ gdy iepªo dopªywa),

T ˆ

^{- znana}temperaturana z± ibrzegu,

ˆ

q

^{- znany strumie« na z± i brzegu,}

q ˆ n

^{- skªadowa normalnastrumienia na brzegu} (narysowana w jednympunk iealejestdanana aªej z± ibrzegu,naktórejniejestznanatemperatura, mawarto±¢

dodatni¡ gdy iepªo dopªywa),

n

^{-wersor normalny nabrzegu.}

przeksztaª i¢ powy»sze sformuªowanie do posta i sªabej potrzebnej w MES nale»y pomno»y¢ obie

strony równania (40)

i

^{przez funk j testow¡ i wy aªkowa¢ przez z± i. Otrzymamy w ten sposób}

odpowiednik zasady pra wirtualny h: znale¹¢ T(x,y), takie »e

T = ˆ T

^na

∂Ω T

ⁱ

Z

Ω

k∇v · ∇T dΩ = Z

Ω

vQ dΩ + Z

∂Ω q

v ˆ q n ds ∀v : v = 0

^na

∂Ω T

⁽⁴¹⁾

Temperatura znana na z± i brzegu jest odpowiednikiem warunku kinematy znego i ze wzgldu

nato,»e wystpuje poza zasad¡ waria yjn¡ nazywa si warunkiempodstawowym(albo Diri hleta),

natomiastznanaskªadowa normalnastrumieniajest wbudowana wfunk jonaª inosi nazw warunku

naturalnego (Neumanna). Tak jak to byªo w zadaniu prta mo»na sformuªowanie waria yjne (41)

zast¡pi¢zadaniem minimaliza jipewnego funk jonaªu.

Rozwa»myzadanieprzykªadowe,którymjestprzepªyw iepªawpryzmaty znympr iezustalon¡

staª¡w zasie iprzestrzeni temperatur¡najednej z± ipobo zni yoraz staªym strumieniem iepªa

napozostaªej z± i pobo zni y. Mo»na zatem rozwa»a¢ tylkojeden przekrój prta. Przyjmijmy,»e

jestto trapez przedstawiony narys. 10z przedstawionymi tam warunkamibrzegowymii ¹ródªem.

PSfrag repla ements

T = 500K ˆ

Q = 120W/m ³

ˆ

q (||ˆ q|| = 100W/m ² ) 1m

1m

2m

Figure10: Przykªadowezadanie2D.Sta jonarnyprzepªyw iepªa. Wspóª zynnikprzewodzenia iepªa

k = 1 _{m deg} ^W

^.

Zastosowanie MES do aproksyma ji rozwi¡zania powy»szego zadania mo»na uzyska¢ w poni»ej

opisany h, poznany h naprzykªadzie zadania1D, etapa h.

1. Sformuªowanie zadania

Skorzystamy ze sformuªowaniawaria yjnego (41). Wiemy,napodstawie analizyfunk jonalnej,

»e jest todobrze postawione zadanie, maj¡ e jednozna zne rozwi¡zanie.

2. Dyskretyza ja - genera ja siatki

Elementy sko« zone w2D mog¡by¢ trójk¡tamialbo zworok¡tami,wogólno± iokrzywolinio-

wy h boka h. Przyjty w przykªadowym zadaniu obszar w ksztaª ie trapezu podzielmy na 3

przystaj¡ e trójk¡typrostok¡tne, równoramienne, bd¡ eelementami sko« zonymi. Generator

siatkidostar za informa jio wspóªrzdny h wzªów:

(0, 1), (0, 0), (1, 1), (1, 0), (2, 0)

^{i o tym}

jakie wzªy deniuj¡ posz zególne elementy sko« zone. Istotna jest kolejno±¢ wzªów przy-

porz¡dkowany h elementom. Zaªó»mynastpuj¡ e rela je: element1 mawzªy(1,2,4),element

2 ma wzªy (4,3,1) i element 3 ma wzªy (3,4,5). Wprowadzony podziaª zyli dyskretyza ja

umo»liwia zdeniowanie funk ji ksztaªtu. S¡ to funk je i¡gªe (dziki temu MES jest sta-

bilna zyli mamy gwaran j, »e zwikszaj¡ li zb elementów sko« zony h zmniejszamy bª¡d

aproksyma ji),przyporz¡dkowane wzªomwtakisposób, »eka»dazfunk jiksztaªtuprzyjmuje

warto±¢1wwybranymw¹leoraz 0wpozostaªy hwzªa h. Wrozwa»anymprzykªadzieka»dy

elementma 3 stopnieswobody a globalnai h li zba (LSS) jest równa 5. Na rys. 12 pokazano

wykres funk ji ksztaªtu zwi¡zanejz wzªem nr. 1.

3. Agrega ja

Stosuj¡ lokaln¡ numera j, wspóª zynniki ma ierzy (

K ^e

^{) i wektora (}

f ^e

^{) elementu}

e

^{obli za}

size wzorów

K _ij ^e = Z

e

k∇ϕ i · ∇ϕ j dΩ, i, j = 1, 2, 3

⁽⁴²⁾

1

2

3

PSfrag repla ements 1

2

3

4 5

0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2

0 0.2 0.4 0.6 0.8 1

500 550 600 650 700

PSfrag repla ements

1

2

3

4

5

Figure 11: Przykªadowe zadanie 2D. Dyskretyza ja i pole temperatury obli zone za pomo ¡ MES

dla3 elementówsko« zony h.

0 1

2 0.5

0.5 1.5

1 1

0.5

0 0 ⁰

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

Figure12: Przykªadowezadanie2D.Wykresfunk jiksztaªtuzwi¡zanejzwzªemnr.1(1/4piramidy).

f _i ^e = Z

e

ϕ i Q dΩ + Z

∂e∩∂Ω q

ϕ i q n ds i, j = 1, 2, 3

⁽⁴³⁾

gdzie

∂e ∩ ∂Ω _q

^{jest z± i¡brzeguelementule»¡ ¡na}"ob i¡»onym"(

q _n

^{jest znane)fragmen ie}

brzegu aªego obszaru.

Warto zauwa»y¢, »e nie ma po prawej stronie zasady waria yjnej skªadnika bez aªki, a za-

tem w globalnym wektorze ob i¡»enia nie bdzie skªadnika obli zanego poza agrega j¡ zyli

odpowiednika wektora

W

^{ukªadu równa«} algebrai zny h (27). Wektoryob i¡»enias¡w ogól- no± i wynikami aªkowania nie tylko po wntrzu elementu, jak to byªo w zadaniu 1D, ale

równie»mog¡ wymaga¢ aªkowania po z± i brzegu

∂Ω q

^.

Zajmijmysielementem1(rys. 13). Przyjmuj¡ ukªadwspóªrzdny hww¹le2,ªatwonapisa¢

1

2 3

Figure13: Przykªadowe zadanie2D. Element nr. 1 z lokalnie ponumerowanymi wzªami.

nastpuj¡ e wzory naponumerowane lokalnie funk je ksztaªtu

ϕ ₁ (x, y) = y

ϕ ₂ (x, y) = 1 − x − y ϕ 3 (x, y) = x

(44)

Ogólnie, funk ja przyjmuj¡ a warto±¢ 1 w w¹le o wspóªrzdny h

(x k , y k )

^{i zeruj¡ a si w}

wzªa howspóªrzdny h

(x i , y i ), (x j , y j )

^{mo»eby¢zapisanawzorem}

ϕ k = A k x + B k y + C k

A k x k + B k y k + C k

,

gdzie

A k , B k , C k

^s¡ wspóª zynnikami równania prostej (

A k x + B k y + C k = 0

⁾ prze hodz¡ ej przez punkty

(x _i , y _i ), (x _j , y _j )

^{. Gradienty funk ji ksztaªtu,} wystpuj¡ e we wzora h na wspóª zynniki ma ierzy sztywno± i, s¡nastpuj¡ ymi wektorami niezale»nymi od

x, y

∇ϕ 1 = [0, 1]

∇ϕ 2 = [−1, −1]

∇ϕ 3 = [1, 0]

(45)

Zatemªatwo obli zy¢, »e dla

k = const K ₁₁ ^e = k R

e [0, 1] · [0, 1] dx dy = k R

e (0 · 0 + 1 · 1) dx dy = kA K ₁₂ ^e = k R

e [0, 1] · [−1, −1] dx dy = −kA K ₁₃ ^e = k R

e [0, 1] · [1, 0] dx dy = 0 K ₂₂ ^e = k R

e [−1, −1] · [−1, −1] dx dy = 2kA . . .

(46)

gdzie

A

^{jest polem elementu} sko« zonego. Ma ierze wszystki h elementów przy zastosowa- nej dyskretyza ji (trójk¡ty przystaj¡ e, z wzªami ponumerowanymi tak aby wzeª przy k¡ ie

prostym byª drugi) s¡jednakowe i równe 9

K ^e = kA





1 −1 0

−1 2 −1 0 −1 1





⁽⁴⁷⁾

9

Elementy sko« zonemog¡ by¢dowolnymi trójk¡tami albo zworok¡tami. Šatwodla ni h uogólni¢ stosowan¡w

tymprostym obli zeniowoprzykªadziemetodobli zaniawspóª zynnikówma ierzyiwektorówelementowy h.

poelemen ie, w nastpuj¡ ysposób.

•

Przyjmujemy zerow¡globaln¡ma ierzsztywno± i

K _nxn

^{izerowyglobalnywektorob i¡»e-}

nia

f _nx1

^{. Poniewa»}

n = 5

K =













0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0













, f =











 0 0 0 0 0













(48)

•

^{Kolejno, dlaka»degoelementuwt.zw. pro edurze} elementowejobli zasima ierz(w tym przykªadzie ju» obli zona)i wektor elementu. I hwspóª zynniki dodajemy doma ierzy i

wektoraglobalnegowmiejs awynikaj¡ eznumerówglobalny hstopniswobodyelementu.

Numery te s¡ adresami (okre±laj¡ kolumny i wiersze) w ma ierzy i wektorze globalnym

podktóre nale»y "dostar zy¢" posz zególne wspóª zynniki ma ierzy i wektora elementu.

(a) Wspóª zynnikiwektoraelementu(wzór(43))s¡wogólno± isumami aªekpowntrzu

i z± i jego brzegu. Obli zmy kolejno potrzebne aªki pamitaj¡ , »e

Q = 120 = const, q _n = 100 = const

¹⁰

Z

e

ϕ 1 Q dx dy = Z

e

ϕ 2 Q dx dy = Z

e

ϕ 3 Q dx dy = 1

3 AQ = 20 Z 1

0

ϕ 1 q n dy = 1

2 q n = 50 Z 1

0

ϕ 2 q n dy = 1

2 q n = 50 Z 1

0

ϕ 3 q n dy = 0

(49)

a zatemwektor ob i¡»eniaelementu 1jest równy

f ⁽¹⁾ =





20 + 50 20 + 50 20 + 0



 =



 70 70 20





⁽⁵⁰⁾

Element 1 mawzªy (1,2,4). Po agrega ji jegoma ierzy i wektora

K = 1 2













1 −1 0 0 0

−1 2 0 −1 0

0 0 0 0 0

0 −1 0 1 0

0 0 0 0 0













, f =











 70 70 0 20

0













(51)

10

R

e ϕ 1 dx dy

^{jestrównaobjto± i ostrosªupajaki tworzywykres funk ji}

ϕ 1

^{(rys. 12), a}

Z ¹

0

ϕ 1 dx

^{jestrównapolu}

trójk¡tajakitworzywykrestejfunk jinaod inkupomidzywzªami 1i2.

f ⁽²⁾ =





20 + 0 20 + 50 20 + 50



 =



 20 70 70





⁽⁵²⁾

Element 2 mawzªy (4,3,1). Po agrega ji jegoma ierzy i wektora

K = 1 2













2 −1 −1 0 0

−1 2 0 −1 0

−1 0 2 −1 0

0 −1 −1 2 0

0 0 0 0 0













, f =











 140

70 70 40 0













(53)

( ) Wektor ob i¡»eniaelementu 3jest równy

f ⁽³⁾ =





20 + 71 20 + 0 20 + 71



 =



 91 20 91





⁽⁵⁴⁾

Element 3 mawzªy (3,4,5). Po agrega ji jegoma ierzy i wektora

K = 1 2













2 −1 −1 0 0

−1 2 0 −1 0

−1 0 3 −2 0

0 −1 −2 4 −1

0 0 0 −1 1













, f =











 140

70 161

60 91













(55)

Otrzymali±myw ten sposób ukªad równa« algebrai zny h

1 2













2 −1 −1 0 0

−1 2 0 −1 0

−1 0 3 −2 0

0 −1 −2 4 −1

0 0 0 −1 1























 T 1

T 2

T 3

T 4

T 5













=











 140

70 161

60 91













(56)

Przyjmuj¡ ozna zenia











 T ₁ T 2

T 3

T ₄ T 5













= u,











 140

70 161

60 91













= f

⁽⁵⁷⁾

zyliwektorystopniswobodyiob i¡»eniamo»nazapisa¢ukªadrówna«algebrai zny hwzwarty

sposób

Ku = f

⁽⁵⁸⁾

Poniewa»

T 2 = T 4 = T 5 = 500

^{, mno»ymy} odpowiadaj¡ e tym stopniom swobody kolumny ma ierzysztywno± iprzezznanewarto± iiprzenosimynapraw¡stron. Nastpnierezygnujemy

zrówna« 2,4,5 zastpuj¡ je równaniami

T 2 = 500, T 4 = 500, T 5 = 500

^{o prowadzi doukªadu}













1 0 −0.5 0 0

0 1 0 0 0

−0.5 0 1.5 0 0

0 0 0 1 0

0 0 0 0 1























 T 1

T 2

T 3

T ₄ T 5













=











 390 500 661 500 500













(59)

któregorozwi¡zaniem s¡

T 1 = 732, T 2 = 500, T 3 = 685, T 4 = 500, T 5 = 500 [K]

⁽⁶⁰⁾

5. Wyniki obli ze«

Znaj¡ stopnieswobody,które przyzastosowany hfunk ja hksztaªtus¡warto± iamirozwi¡za-

nia przybli»onego w wzªa h mo»na opra owa¢ wyniki element po elemen ie. Wiedz¡ , na

podstawie rela jielement-wzªy jakie s¡ numery wzªów pobieramy wspóªrzdne ty h wzªów

orazstopnie swobody elementu (ozna zaj¡ je zzastosowaniemlokalnejnumera ji

T ₁ , T ₂ , T ₃

⁾ⁱ

zapisujemy rozwi¡zanieMES w rozwa»anym elemen ie sko« zonym

T _h ^e (x) = T 1 ϕ 1 (x, y) + T 2 ϕ 2 (x, y) + T 3 ϕ 3 (x, y)

⁽⁶¹⁾

Przykªadowo, dla elementu 1, stosuj¡ numera j lokaln¡

T 1 = 732, T 2 = 0.500, T 3 = 500

^{, a}

zatem

T _h ⁽¹⁾ (x, y) = 732y + 500(1 − x − y) + 500x = 232y + 500

⁽⁶²⁾

Stosuj¡ prawo Fourieramo»emy obli zy¢strumie« iepªa wtym elemen ie

q ⁽¹⁾ _h = −k∇T _h ⁽¹⁾ = [0, −232]

⁽⁶³⁾

Wektory

q _h

^{s¡ staªe w elementa h i nie i¡gªe midzy nimi. Podobnie mo»na obli zy¢}

q ⁽²⁾ _h = [47, −185]

ⁱ

q ⁽³⁾ _h = [0, −185]

^{. Wykres i¡gªego pola} temperatury pokazano na rys. 11b, a wektory strumienia iepªa obli zonego na podstawie przyjtej bardzo zgrubnej dyskretyza ji

narys. 14.

6. Wiarygodno±¢ wyników

Zastosowana dyskretyza ja byªa bardzo zgrubna. Aby osza owa¢ bª¡d uzyskany h wyników

nale»aªoby powtórzy¢ obli zenia z przynajmniej dwa razy mniejszymi elementami o daªoby

ok. 4x wi ej stopniswobody. Porównanietakuzyskanegowyniku z przedstawionym powy»ej

daªoby informa j o bªdzie tego pierwszego rozwi¡zania. Na rys. 15 porównano wyniki

uzyskane przy 5 i 120 stopnia h swobody. Poniewa» rozwi¡zanie jest do±¢ wolno zmieniaj¡ ¡

sifunk j¡ oba wynikinieró»ni¡ sizna znie.

0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 -0.2

0 0.2 0.4 0.6 0.8 1 1.2

Figure14: Przykªadowe zadanie 2D. Wektoryintensywno± i strumienia iepªa obli zone zapomo ¡

MESdla3 elementów sko« zony h.

Na zako« zenie tego rozdziaªu warto wspomnie¢, »e stosowane tu równanie Lapla e'a (± i±lej Pois-

son'a) modeluje równie» inne zjawiska ni» sta jonarny przepªyw iepªa. Przykªadowo, zamiast tem-

peratury pozwala obli zy¢ ustalone st»enie substan ji w zjawisku dyfuzji, poten jaª elektry zny

wywoªany ªadunkiem elektrostaty znym, ugi ia membrany (na i¡gnitego pªótna) od ob i¡»enia

prostopadªego do jej powierz hni, deforma j przekroju prta pryzmaty znego poddanego skr a-

niu.... Dla wszystki h ty h zada« MES pozostaje bez zmian, gdy» wszystkie wspomniane

zjawiska s¡opisywane takimsamym równaniem.

8 MES dla liniowej spr»ysto± i

Naj z± iej rozwi¡zywanym zadaniem przy projektowaniu konstruk ji jest statyka iaª spr»ysty h.

Dlatego zajmijmysi wybranymi sz zegóªamizastosowania MESdo tegotypuzada«. Za znijmy od

sformuªowania matematy znego znanego z me haniki pozwalaj¡ ego obli zy¢ niezmienne w zasie

przemiesz zenia

u (x)

^, odksztaª enia

ǫ (x)

^{i napr»enia}

σ (x)

^{w zakresie liniowym, przy zerowym}

przyspieszeniu, gdzie

x

^{ozna za wektor} wspóªrzdny h przestrzenny h punktu



 

 



 

 



−divσ = q w Ω

σ = Cǫ w Ω

ǫ = sym(∇u) w Ω u = 0 na ∂Ω u

σn = ˆt na ∂Ω t

(64)

q, ˆt

^{s¡ znanymi} ob i¡»eniami wewn¡trz i na brzegu,

C

^{ozna za tensor parametrów} materiaªowy h reprezentowany przez 4 wska¹nikow¡ ma ierz(zwan¡ zasem skrótowo ma ierz¡sztywno± i),

n

^jest

skierowanym na zewn¡trz wersorem prostopadªym do brzegu (rys. 16). Sformuªowaniem sªabym

zadaniaspr»ysto± i mo»e by¢ zasada pra wirtualny h: znale¹¢

u

⁽

x

^{), takie »e}

u = 0

^na

∂Ω u

ⁱ

Z

Ω

ǫ (v) : Cǫ(u) dΩ = Z

Ω

v · q dΩ + Z

∂Ω t

v · ˆt ds ∀v : v = 0

^na

∂Ω u

⁽⁶⁵⁾

0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 0

0.2 0.4 0.6 0.8 1

500 550 600 650 700

Figure 15: Przykªadowe zadanie 2D. Rozwi¡zania MES uzyskane przy 5 i 120 stopnia h swobody.

Maksymalne temperatury ró»ni¡ si o mniej ni» 1 % ale rozkªady temperatur wewn¡trz obszaru

wykazuj¡ wiksze ró»ni e.

PSfrag repla ements

∂Ω _t

∂Ω _u

q

ˆt n

Figure16: Ozna zeniastosowane wzadaniuspr»ysto± i,2D.

q

^{-ob i¡»eniewewn¡trz iaªa,}

ˆt

^-znane

ob i¡»eniepowierz hniowena brzegu,

n

^{-wersor normalny nabrzegu.}

Problem(65)mo»eby¢sformuªowanyrównie»jakominimumfunk jonaªu,t.zn. znale¹¢

u(x), u = 0

^na

∂Ω u

^{, takie »e}

J(u) = min v : v =0 na ∂Ω u

J(v)

⁽⁶⁶⁾

gdzie

J

^{jest aªkowit¡ energi¡}poten jaln¡ukªadu

J(v) = 1 2

Z

Ω

ǫ(v) : Cǫ(v) dΩ − Z

Ω

v · q dΩ − Z

∂Ω t

v · ˆt ds

⁽⁶⁷⁾

Podstawow¡ e h¡zadaniame hanikijestfakt»eniewiadomafunk japrzyjmujewarto± iwektorowe.

Zastanówmy si jak mo»e wygl¡da¢ aproksyma ja rozwi¡zania którym s¡ tutaj przemiesz zenia

maj¡ e w ogólno± i trzy skªadowe

u x , u y , u z

^{. Dla} uprosz zenia zajmijmy si pªaskim stanem od- ksztaª enia, bd¡ ego dobrym modelem dla obiektów pryzmaty zny h z ob i¡»eniem i podpar iem

niezale»nym od poªo»enia w kierunku osi obiektu (np. mur oporowy) 11

. Ponownie zajmijmy si

11

Pªaski stannapr»eniaanalizujesi wbardzo podobny sposób i mazastosowaniedo obiektówomaªym jednym

wymiarze(grubo± i),np. ± ian.