23
Uvod do Metody Konecnych Prvku (MKP)
I Vetsina fyzikalnıch jevu muze byt popsana systememparcialnıch diferencialnıch rovnic.
I Analyticke resenı klasickymi metodami pro obecne oblastije velmi obtızne ci prakticky nemozne.
I MKP (Finite Element Method - FEM1) je nejcasteji uzıvana,systematicka a univerzalnı metoda pro numericke resenıproblemu.
1Google “FEM” > 14 mil. odkazu, > 500 knih o MKP, vydaje > 1 mld.USD na FEM Software a vypocetnı cas
Idea MKP
Historie MKP1943 Courant - aplikace variacnıch principu, polozil zaklad
matematicke teorie MKP1950 prvnı inzenyrska aplikace v letectvı, Boeing&Bell.
M.J.Turner, R.W.Clough (→ Berkeley), M.C.Martinpublikovali jeden z prvnıch clanku.
I Berkeley: E. Wilson, R.L.Taylor a jejich PhD studenti:T.J.R. Hughes, C. Felipa, K.J. Bathe
I Swansea: O.C. Zienkiewicz, B. Irons, R.Owen1960 E. Wilson - prvnı MKP program (freeware)1965 Nastran; 1969 Ansys (hodnota spolecnosti ≈ 1.8 mld. $);
1978 Abaqusdnes MKP aplikovana pro resenı soucasnych
vedecko-technickych problemu - komplexnı navrh letadel,simulace vyrobnıch procesu, narazove testy automobilu,navrh spalovacıch motoru, ochranych obalek jadernychreaktoru, seismicka analyza prehrad, aplikace v lekarstvı -simulace provadenı plastickych operacı, ...
Historicky vyvoj rychlosti pocıtacu
Z obrazku je patrny linearnı narust rychlosti (logaritmicke merıtko). Toho si prvne vsiml G. Moore2 v roce 1965 a
formuloval empiricke pravidlo, zname jako Mooruv zakon:Slozitost soucastek (pocet tranzistoru na cipu a jejich
vykon) se kazde dva roky zdvojnasobı pri zachovanı stejne ceny.
2Spoluzakladatel firmy Intel (1968)
Historicky vyvoj vybranych cen v USA3
1968 2005CDC 6600 (0.5-1 Mflopf) $8000000512 Beowulf cluster, 2003 (1 Tflop) $500000PC (200-1600 Mflops) $500-$3000MSc Engineer, starting salary $9000 $51000Assistant Prof. $11000 $55000Tuition at Northwestern (1 year) $1800 $32000GM, Ford sedan $3000 $22000
3bez vlivu inflace; prevzato z Hughes-Belytschko FEM course
Top500 Supercomputers (www.top500.org)
Poradı v cervnu 2019
Rank Site System Cores Rmax4 Power(kW)
1 DOE/SC/Oak RidgeNational Lab
Summit - IBM Power System AC922,IBM POWER9 22C 3.07GHz, NVIDIAVolta GV100, Dual-rail Infiniband
2,414,592 148,600 10,096
2 DOE/NNSA/LLNL Sierra-IBM Power SystemS922LC,IBM POWER9 22C 3.1GHz,NVIDIA Volta GV100, Infiniband
1,572,480 95,640 7,438
3 National Supercom-puting Center inWuxi, China
Sunway TaihuLight - Sunway MPP,Sunway SW26010 260C 1.45GHz
10,649,600 93,014 15,371
282 IT4Innovations VSB-Ostrava
Salomon - SGI ICE X, Xeon E5-2680v3 12C 2.5GHz, Infiniband, Cen-tOS
76,896 1,457.73 4,806
4Maximalnı vykon v testu LINPACK (TFlops/s)
Prıklady Aplikace MKP
MKP: Prıme odvozenı pro tazene-tlacene pruty
I Postup analogicky odvozenı deformacnı metodyI Spocıva ve vyjadrenı uzlovych sil v zavislosti na koncovych
posunech prvkue
x
F , u1 1
eeF ,u2 2
e e
1 2
Mame k dispozici nasledujıcı rovnice:I Podmınky rovnovahy mezi vnitrnımy silami (σ) a uzlovymi
silami (F e1 ,F
e2 ): F e
1 = −σA, F e2 = σA
I Vztah mezi napetım a deformacı (Hookuv zakon) σ = EεI Definice deformace jako pomerneho protazenı prutuε = ∆l
l
Pro napetı tedy platı
σ = Eε = E∆ll
= Eue
2 − ue1
l
A pro koncove sıly konecne dostavame
F e1 = −σA =
EAl
(ue1 − ue
2); F e2 = σA =
EAl
(ue2 − ue
1)
Maticove to lze zapsat nasledovne:{F e
1F e
2
}=
[k −k−k k
]{ue
1ue
2
}, F e = K ere
kde k =EAl
.
{F e
1F e
2
}=
[k −k−k k
]{ue
1ue
2
}I Pokud ue
1 = ue2 (posunutı prutu jako tuheho celku)
nevzniknou zadne vnitrnı sıly a F e1 = F e
2 = 0I Matice tuhosti je symetricka a singularnıI Linearita je dusledkem predpokladu (linearita
konstitutivnıch vztahu, geomerickych a rovnovaznychrovnic)
Sestavenı matice tuhosti konstrukce - Lokalizace
������������������������������������
������������������������������������
l
E , A
E , A
1 2
1 1
2 2
l
1 2f = 10 f = 5
f ,u 1 2 3
r ,u3
f ,u 21
1 3(1) (2)2
Pro kazdy prut mame kdispozici vztah mezikoncovymi silami a posuny
F 1 = K 1r1, F 2 = K 2r2
Podmınky rovnovahy vuzlech:
(1)←: F 11 − f1 = 0
(2)←: F 12 + F 2
1 − f2 = 0(3)←: F 2
2 − r3 = 0
F12
F12 F
22F1
1f1 f2 r3
F22F1
2F12F1
1
11u
21u
12u
22u
f ,u 1 2 3
r ,u3
f ,u 21
1 3(1) (2)2
1 2 3
u11 = u1,u1
2 = u21 = u2,u2
2 = u3
Uvazıme-li kompatibilitu posunutı mezi prvky, muzeme vyjadritkoncove sıly pomocı uzlovych posunu:
(1)←: (k1u1 − k1u2)− f1 = 0(2)←: (−k1u1 + k1u2) + (k2u2 − k2u3)− f2 = 0(3)←: (−k2u2 + k2u3)− r3 = 0
Maticove zapsano:
(1)←:(2)←:(3)←:
F 11
F 12
0
+
0F 2
1F 2
2
=
f1f2r3
Abychom mohli koncove sıly vyjadrene na jednotlivych prutechsnadno scıtat, vyjadrıme je pomocı globalnıch posunutı{
F 11
F 12
}= K e
1
{ue
1ue
2
}= K e
1
{u1u2
}{
F 21
F 22
}= K e
2
{ue
1ue
2
}= K e
2
{u2u3
}
Lokalnı vektory koncovych sil a posunutı a matici tuhostirozsırıme pridanım nulovych prvku tak, aby obsahovali vsechnyhodnoty pro celou konstukci:
F 11
F 12
0
︸ ︷︷ ︸˜F
1
=
k1 −k1 0−k1 k1 0
0 0 0
︸ ︷︷ ︸
˜K1
u1u2u3
0F 2
1F 2
2
︸ ︷︷ ︸˜F
2
=
0 0 00 k2 −k2
0 −k2 k2
︸ ︷︷ ︸
˜K2
u1u2u3
Pak muzeme podmınky rovnovahy elegantne zapsat jako
F1
+ F2
= f
Po dosazenı za koncove sıly pak k1 −k1 0−k1 k1 0
0 0 0
︸ ︷︷ ︸
˜K1
u1u2u3
︸ ︷︷ ︸r
+
0 0 00 k2 −k2
0 −k2 k2
︸ ︷︷ ︸
˜K2
u1u2u3
︸ ︷︷ ︸r
=
f1f2r3
︸ ︷︷ ︸f
Maticove pak(K
1+ K
2)︸ ︷︷ ︸
K
r = f
Tento proces lze formalizovat. Definujme pro kazdy prvek tzv.distribucnı matici Le, tak ze platı re = Ler . Napr. pro prvek 1platı:
r1 =
{u1
1u1
2
}=
[1 0 00 1 0
]︸ ︷︷ ︸
L1
u1u2u3
= L1r
To platı obdobne i pro vektor prave strany:
F 1 = L1F1
Nase rovnice na prvku K ere = F e pak muzeme vyjadrit vglobalnıch slozkach jako
K eLer = F e
A prenasobenım zleva maticı LeT nakonec dostavame
LeT K eLe︸ ︷︷ ︸˜K
e
r = LeT F e︸ ︷︷ ︸f
Lokalizace - algoritmizace sestavenı globalnı maticetuhosti
I Zavedenı kodovychcısel, reprezentujıcıocıslovanı neznamychposunu r = {u1,u2,u3}T
I Jejich prirazenı lokalnımposunum a silam naprvcıch
I Prvky lokalnı maticetuhosti se pricıtajı doprvku globalnı maticeurceneho kodovymicısly.
+
1
2
3
1
2
2
3
2 31
1 2 2 3
1
2
3 k
2 k
2−
k1
k1
−
k1
−1
k k2
k2
−
1 k
k1
k1
−
k1
− u2
1
u1
1
k2
k2
− u1
2
u2
2
k2
k2
−
1
2
3u
u
u
=
=
1
=
F11f1
F12
F12
f2
F22
r3
F22F1
2F12F1
1
F12
F12
f2F11f1 F
22
r3
F11
F12
F12
F22
1
2
3
f
f
r
1
2
3
f ,u 1 2 3
r ,u3
f ,u 21 (1) (2)
1 2 3
2 3
Lokalizace - algoritmizace sestavenı globalnı maticetuhosti
I Zavedenı kodovychcısel, reprezentujıcıocıslovanı neznamychposunu r = {u1,u2,u3}T
I Jejich prirazenı lokalnımposunum a silam naprvcıch
I Prvky lokalnı maticetuhosti se pricıtajı doprvku globalnı maticeurceneho kodovymicısly.
+
1
2
3
1
2
2
3
2 31
1 2 2 3
1
2
3 k
2 k
2−
k1
k1
−
k1
−1
k k2
k2
−
1 k
k1
k1
−
k1
− u2
1
u1
1
k2
k2
− u1
2
u2
2
k2
k2
−
1
2
3u
u
u
=
=
1
=
F11f1
F12
F12
f2
F22
r3
F22F1
2F12F1
1
F12
F12
f2F11f1 F
22
r3
F11
F12
F12
F22
1
2
3
f
f
r
1
2
3
f ,u 1 2 3
r ,u3
f ,u 21 (1) (2)
1 2 3
2 3
Lokalizace - algoritmizace sestavenı globalnı maticetuhosti
I Zavedenı kodovychcısel, reprezentujıcıocıslovanı neznamychposunu r = {u1,u2,u3}T
I Jejich prirazenı lokalnımposunum a silam naprvcıch
I Prvky lokalnı maticetuhosti se pricıtajı doprvku globalnı maticeurceneho kodovymicısly. +
1
2
3
1
2
2
3
2 31
1 2 2 3
1
2
3 k
2 k
2−
k1
k1
−
k1
−1
k k2
k2
−
1 k
k1
k1
−
k1
− u2
1
u1
1
k2
k2
− u1
2
u2
2
k2
k2
−
1
2
3u
u
u
=
=
1
=
F11f1
F12
F12
f2
F22
r3
F22F1
2F12F1
1
F12
F12
f2F11f1 F
22
r3
F11
F12
F12
F22
1
2
3
f
f
r
1
2
3
f ,u 1 2 3
r ,u3
f ,u 21 (1) (2)
1 2 3
2 3
Kodova cısla tedy predstavujı mapovanı mezi lokalnımcıslovanım na prvku a globalnım cıslovanı pro celou konstukciFormalne pak pıseme K g =
∑i K
i
Obdobny proces i pro vektory
Zavedenı okrajovych podmınek
I Matice tuhosti konstrukce je singularnı - obsahuje posunutıtelesa jako tuheho celku
I Zavedenım okrajovych podmınek dojde k regularizacisoustavy
I Rozdelme podmınky rovnovahy na ty, ktere prıslusı volnyma predepsanym stupnum volnosti[
K uu K upK up K pp
]{uu
}=
{fr
}I Z prvnıho radku muzeme spocıst neznamy vektor posunu:
u = K−1uu (f − K upu)
I A z druhe rovnice muzeme nasledne dopocıtat neznamereakce: r = K upu + K ppu.
