Institutionen för naturvetenskap och teknik
En matematisk modell av
släggkastets rotationsfas
Caroline Hermansson
Anton Karlsson
Örebro universitet
Institutionen för naturvetenskap och teknik Matematik C
En matematisk modell av släggkastets
rotationsfas
Caroline Hermansson Anton Karlsson
Juni 2014
Handledare: Mårten Gulliksson och Henrik Andreasson Examinator: Marcus Sundhäll
Självständigt arbete, 15 hp Matematik, Cnivå
Sammanfattning
Med vår uppsats visar vi hur det med hjälp av Lagranges metod går att ta fram en matematisk modell av rotationsfasen i ett släggkast. Släggkastaren och släggan utgörs i uppsatsen av en sammansatt stelkroppsmodell. Vi tar fram kastarens rörelseekvationer vilka vi får genom Lagranges ekvation. I uppsatsen har vi härlett Lagranges ekvation genom Newtons andra lag och Hamiltons princip. Resultatet av modellen presenteras i form av en plot i Matlab av slägghuvudets position under olika tider under rörelsen.
Innehåll
1 Introduktion 5 1.1 Inledning . . . 5 1.2 Bakgrund . . . 5 1.2.1 Kinetisk energi . . . 6 1.2.2 Potentiell energi . . . 9 1.2.3 Eulervinklar . . . 102 Lagranges funktion och ekvation 13 2.1 Härledning av Lagranges ekvation genom Newtons andra lag . 14 2.1.1 Satsen om generellt koordinatbyte . . . 17
2.2 Hamiltons princip . . . 19
2.3 Ett exempel på Lagranges metod . . . 21
3 Släggkastets rotationsfas utifrån Lagranges metod 26 3.1 Beräkning av rörelseekvationer för vår modell . . . 26
3.1.1 Cylinderns kinetiska energi . . . 28
3.1.2 Punktmassans kinetiska energi . . . 29
3.1.3 Cylindern och punktmassans potentiella energi . . . . 30
3.1.4 Lagrangefunktionen för kroppen . . . 30
3.1.5 Rörelseekvationer . . . 30
3.2 Lagranges ekvationer på matrisform . . . 31
3.3 Implementering i Matlab av matrisformuleringen . . . 32
3.4 Avancerad simulering . . . 33
3.5 Resultat . . . 33
4 Diskussion 38
A Matlabkod för ett exempel på Lagrange metod 42
B Matlabkod för släggkastets rotationsfas 44
Kapitel 1
Introduktion
1.1 Inledning
Tanken med uppsatsen är att bygga upp en matematisk modell av rotations-fasen som utförs i friidrottsgrenen slägga, utifrån en sammansatt stelkropps-modell. Genomgående för uppsatsen är att vår modell samt alla exempel kommer att utgå ifrån att kroppen be nner sig i R3 om inte annat anges.
Läsaren antas vara matematikstuderande på universitet- eller högskolenivå. Uppsatsen kan även vara av intresse för den friidrottsintresserade som vill se hur den matematiska teorin bakom släggkastets rotationer kan se ut.
Vi har valt att strukturera uppsatsen på så sätt att vi delat in den i fyra kapitel. Första kapitlet innehåller en bakgrund där de viktigaste begreppen de nieras och där vi använder exempel för de de nitioner vi anser kräver ytterligare beskrivning. I kapitel två beskrivs vilken metod vi använder oss av och varför vi använder den, samt visar hur metoden kan härledas från Newtons andra lag och genom Hamiltons princip. Tredje kapitlet innehåller vår stelkroppsmodell där vi med hjälp av Lagranges metod tar fram rörelse-ekvationer för vår valda kropp. Rörelserörelse-ekvationerna är kopplade, icke-linjära andra ordningens ordinära dierentialekvationer. I detta kapitel nns även resultatet av vår modell där slägghuvudets rörelse under rotationsfasen i ett släggkast visualiseras i form av en plot från Matlab. Kapitlet tar även upp mer avancerad modellering och simulering. Det sista kapitlet innehåller en diskussion där vi bland annat jämför vår uppsats med en tidigare studie inom samma område.
1.2 Bakgrund
Inom den klassiska mekaniken, som grundar sig i Isaac Newtons rörelselagar, åter nns den analytiska mekaniken. Denna inrikning kom under 1700-talet främst att utvecklas av Leonhard Euler och Joseph-Louis Lagrange och det är framförallt i deras teorier vi grundar vår uppsats. Tröghetsmoments- och
dierentialkalkylsberäkning är exempel på begrepp som används inom den analytiska mekaniken och som vi kommer att använda [2].
Inom idrotten kan den analytiska mekaniken, som en del av biomekanik, komma till hjälp vid rörelseanalys och utveckling av teknik i speci ka gre-nar [5]. Det är här vi tar vår utgångspunkt och tanken är att ta fram en matematiskt modell av rotationsfasen i ett släggkast.
För att få en bild av hur släggkastet ser ut följer nedan en kort besrkivning av dess rörelser;
Släggkastaren står längst bak i kastringen med ryggen mot kastrikt-ningen och håller i släggan så att den är placerad på marken framför kroppen.
För att påbörja kastet och få upp hastigheten på släggan gör kastaren ett par försvängar, det vill säga att armarna roteras i axelhöjd så att släggan i sin tur roteras kring kroppen.
Därefter påbörjas rotationerna i fötterna och benen så att kastaren ro-terar 4.5 varv, samt för yttas i en rak linje bakåt mot utkastriktningen, innan släggan släpps.
Genom att se människokroppen som en tredimensionell stel kropp, det vill säga en kropp som ej är deformerbar [4], fås en modell att utgå ifrån. Denna kommer att rotera likt rotationerna i den tredje ovannämnda punkten av släggkastet. Modellen kommer så småningom, genom Lagrangefunktionen, att utvecklas till Lagranges ekvationer som sammansatta rörelseekvationer. Lösningen till dessa ekvationer beskriver den stela kroppens läge i varje ögon-blick av rörelsen [4].
Vi vill lyfta fram och förtydliga vissa geometriska begrepp som vidare kommer att användas. Vi kommer hädanefter att beteckna det rums xa ko-ordinatsystemet som OXY Z och det lokala koko-ordinatsystemet som Oxyz, se gur 1.1, för en partikel [17].
1.2.1 Kinetisk energi
Som referens till de de nitioner som ges här, hänvisar vi till [4]. Då vi så småningom vill beräkna rörelseekvationerna för vår framtagna modell krävs att tröghetsmomentet beräknas. Denna storhet beskriver ett partikelsystems förmåga att ändra sin rotationshastighet och de nieras för ett partikelsystem som
IQ =X
k
mkr2k
där Q betecknar partikelsystemets momentaxel, det vill säga den axel som partikelsystemet roterar kring. rkär det vinkelräta avståndet från varje
Figur 1.1: Lokalt koordinatsystem, Oxyz, i en partikel i ett globalt koordi-natsystem, OXYZ.
i kroppen ges tröghetsmomentet istället av IQ=
Z
r(x, y, z)2dm
där r(x, y, z) är det vinkelräta avståndet från (x, y, z) till Q och där dm = ρ(x, y, z) dV och ρ(x, y, z) är partikelsystemets densitet och dV är varje litet volymselement.
Tröghetsmomenten för ett partikelsystem kring x-, y- respektive z- axeln i ett tredimesionellt koordinatsystem de nieras som
Ixx =X k mk(y2k+ zk2), Iyy= X k mk(zk2+ x2k), Izz= X k mk(x2k+ y2k).
Om kroppen har en kontinuerlig massfördelning ersätts summorna med
Ixx = Z y2+ z2dm, Iyy = Z z2+ x2dm, Izz = Z x2+ y2dm där dm = ρ(x, y, z) dV = ρ(x, y, z) dx dy dz.
Låt nu en cylinder med jämt fördelad massa utgå från origo i Oxyz och rotera denna kring z - axeln. Cylindern antas ha radien r, höjden l1, massan
m1 och volymen V , se gur 1.2. För denna rotation ges tröghetsmomentet
av Izz = Z x2+ y2dm = Z Z Z x2+ y2 ρ dx dy dz = ρ Z Z Z x2+ y2dx dy dz. Efter beräkningar av integralen, där vi antar konstant densitet ρ(x, y, z) = m/V , får vi
Izz =
l21m1
Figur 1.2: Cylinder som roterar kring Z− axeln. Där r = radien, m = massan, h =höjden och V = volymen.
Beräkningen av Izz kommer att användas som en del i beräkningarna av
tröghetsmomenten för vår modell i kapitel 3.
Eftersom att vår kommande modell utgår ifrån ett tredimesionellt koor-dinatsystem med x-, y- respektive z - axel kommer det även att uppstå så kallade tröghetsprodukter [4]. Dessa ger utöver tröghetsmomenten ytterligga-re information om massfördelningen i kroppen. Om någon tröghetsprodukt är skild från noll betyder det att det är asymmetri i massfördelningen vid de två axlarna som produkten utgörs av.
För rotation kring t.ex. x - respektive y - axeln de neras tröghetspro-dukten för ett partikelsystem som
Ixy = Iyx= −
X
k
mixiyi.
För en kontinuerlig massfördelning i kroppen har vi Ixy = Iyx= −
Z
xy dm.
För cylindern i ovanstående exempel skulle tröghetsprodukten som upp-står kring x- respektive y - axeln de nieras enligt
Ixy = Iyx= − Z xy dm = − Z Z Z xyρ dx dy dz = −ρ Z Z Z xy dx dy dz. Efter beräkningar av integralen får vi
Ixy = Iyx= 0.
Tröghetsmomenten, tillsammans med tröghetsprodukterna, bildar en trög-hetstensor. Denna utgörs av en 3×3 matris som består av tröghetsmomenten
kring x-, y- och z - axeln och tröghetsprodukterna. Matrisen är alltid sym-metrisk och har 6 oberoende element,
I = Ixx Ixy Ixz Iyx Iyy Iyz Izx Izy Izz .
Energin som en kropp har på grund av en viss rörelse kallas för kinetisk energi. Den kinetiska energin utgörs dels av rotationsenergi, dels av transla-tionsenergi. Translationsenergi är den energi som kroppen använder vid en linjär rörelse, det vill säga en rörelse utan rotation [17]. För att kunna ta fram Lagrangefunktionen krävs att den kinetiska energin beräknas. Till den-na krävs vinkelhastigheten, ω = dθ/dt där θ(t) är rotationsvinkeln och där även dθ/dt kan uttryckas som ˙θ [16]. Om rotationerna sker kring er än en axel har vi generellt ω = [ωx, ωy, ωz]T, där ωx, ωy och ωz är
vinkelhastighe-terna kring respektive axel.
Antag att vi har en kropp som rör sig fritt i OXY Z vars tyngdpunkt har translationshastigheten v = [vx, vy, vz]T och rotationen kring
tyngd-punkten är ω. Translationsenergin och rotationsenergin de nieras då som Ttra = mvTv/2 och Trot = ωTIω/2 där m är kroppens massa. Genom att
addera dessa två uttryck fås den sammanlagda formeln för den kinetiska energin enligt T = Ttra+ Trot= 1 2mv Tv +1 2ω TIω. (1.1)
Om vi istället har en punktmassa med massan m blir den kinetiska ener-gin T = 1 2mvp Tv p, (1.2) där vp= v + w × r (1.3)
och r = [x, y, z]Tär positionsvektorn för punktmassan i förhållande till origo.
1.2.2 Potentiell energi
Genomgående för uppsatsen är att alla exempel och beräkningar utgår ifrån att kropparna be nner sig i ett konservativt kraftfält. För att de niera det konservativa kraftfältet antar vi att vi har en partikel (eller kropp) som påverkas av en kraft F (r), r = r(x, y, z). Den potentiella energin i en punkt brelativt en punkt a de nieras i sin tur som arbetet att föra en partikel (eller kropp) från punkt a till punkt b. Om det existerar en potential V sådan att F = ∇V sägs kraftfältet vara konservativt. Följande sats, härledd ifrån [7], gör det möjligt att beräkna potentialen för ett konservativt kraftfält vilket vi senare kommer att behöva.
Sats 1.2.1. Låt F vara ett kraftfält med potentialen V , det vill säga F = ∇V, i det öppna området Ω. För varje kurva γ i Ω gäller då att
Z
γ
F · dr = V (a) − V (b) där a och b är start- respektive slutpunkt för γ.
Bevis. Låt r = r(t), α ≤ t ≤ β, vara en parameterframställning av γ. Speci-ellt är
r(α) = a, r(β) = b. Enligt kedjeregeln är
d
dtV (r(t)) = ∇V (r(t)) · ˙r(t) = F (r(t)) · ˙r(t). Utifrån detta får vi att
Z γ F · dr = Z β α F (r(t)) · ˙r(t) dt = Z β α d dtV (r(t)) dt = V (r(α)) − V (r(β)) = V (a) − V (b).
Notera här att potentialen ej är unik ty en godtycklig konstant C kan adderas till V utan att värdet efter beräkning av ∇V påverkas.
Det är även så att energin i ett konservativt kraftfält bevaras enligt kon-serveringslagen. Detta innebär att T +V är konstant och att systemet vinner lika mycket kinetisk energi T som det förlorar potentiell energi V [6].
För att ta fram den fullständiga Lagrangefunktionen krävs även den po-tentiella energin. Tyngdkraften på en kropp är F = −mgez där ez är
en-hetsvektorn i z - led, m är kroppens massa samt g är gravitationsacceleratio-nen. Genom att beräkna integralen i satsen ovan givet F får vi potentialen V = mgz där z utgörs av höjden över xy - planet. Vi har här antagit att att V (0) = 0. Följdaktligen är F konservativt och på så sätt kan vi skriva den potentiella energin som
V = mgz. (1.4)
1.2.3 Eulervinklar
När vi i kapitel 3 presenterar vår modell kommer positionen för kroppen att beskrivas med hjälp av Eulervinklar. Dessa beskriver förhållandet mellan det lokala koordinatsystemet Oxyz för kroppen och det globala, rums xa, koordinatsystemet OXY Z. För att beskriva rotationen mellan koordinatsy-stemen används tre rotationsvinklar, t.ex. θ1, θ2 och θ3, som i sin tur kallas
X-y'-x, X-z'-x, Y-x'-y, Y-z'-y, Z-x'-z, Z-y'-z, där ' och står för varje ny axel som bildas efter föregående vridning.
Varje rotation beskrivs av en 3 × 3 matris varpå det nns en sådan till varje axel. Dvs en matris om rotationen sker kring en x - axel, en matris om rotationen sker kring en y - axel samt en matris för rotation kring en z - axel [9].
Vi ska nu beskriva de s.k klassiska Eulervinklarna. Välj OXY Z och Oxyz samt i vilken ordning rotationerna ska ske. I detta fall väljer vi rotationer moturs kring Z-x0-z00, se gur 1.3. Rotationen kring Z - axeln de nieras
genom vinkeln θ1, vilket ger oss ett nytt koordinatsystem Ox0y0z0. Rotationen
ser ut enligt följande
Rz(θ1) = cos θ1 sin θ1 0 − sin θ1 cos θ1 0 0 0 1 . (1.5)
Utifrån detta görs rotationen kring x0 - axeln genom vinkeln θ
2, vilket ger
koordinatsystemet Ox00y00z00 och rotationen sker enligt
Rx(θ2) = 1 0 0 0 cos θ2 sin θ2 0 − sin θ2 cos θ2 . (1.6)
Till sist görs rotationen kring z00 - axeln genom vinkeln θ
3, vilket ger oss det
sökta koordinatsystemet Oxyz. Denna rotation görs enligt
Rz(θ3) = cos θ3 sin θ3 0 − sin θ3 cos θ3 0 0 0 1 . (1.7)
En anmärkning i avseende med ovanstående rotationsmatriser är att de gäller vid rotation moturs av ett koordinatsystem xt i en kropp, i ett glo-balt koordinatsystem. En annan möjlighet är att tänka sig en vektor, xt i en kropp som moturs roterar i ett globalt koordinatsystem, vilket i sin tur genererar transponatet av rotationsmatriserna, se [13].
För att beräkna rotationsenergin i den kinetiska energin för en kropp krävs även vinkelhastighetsvektorn ω = [ωx, ωy, ωz]T, där ωx, ωy och ωz är
vinkelhastigheten kring respektive axel. Vid rotation med Eulervinklar gäller enligt [17] att ωx ωy ωz = s2s3 c3 0 s2c3 −s3 0 c2 0 1 ˙ θ1 ˙ θ2 ˙ θ3 = ˙ θ1s2s3+ ˙θ2c3 ˙ θ1c3s2− ˙θ2s3 ˙ θ1c2+ ˙θ3 (1.8)
där s1 = sin(θ1), s2 = sin(θ2), s3 = sin(θ3), c1 = cos(θ1), c2 = cos(θ2) och
Figur 1.3: Eulervinklar. θ1, θ2 och θ3 är vridningar moturs kring Z, x'
re-spektive z - axeln.
Vi kan även ta fram en rotationsmatris genom produkterna av (1.5), (1.6) samt (1.7) vilket ger oss [8]
R = Rz00(θ3)Rx0(θ2)RZ(θ1) = c1c3− c2s1s3 −c1s3− c2c3s1 s1s2 c3s1+ c1c2s3 c1c2c3− s1s3 c1s2 s2s3 c3s2 c2 .
Kapitel 2
Lagranges funktion och
ekvation
Den funktion som inom den klassiska mekaniken används för att förklara rö-relseekvationerna för ett konservativt mekaniskt system kallas för Lagrange-funktionen. Funktionen uttrycks som dierensen mellan systemets kinetiska och dess potentiella energi, det vill säga L = T − V . Genom att sätta in Lagrangefunktionen i Lagranges ekvationer fås systemets rörelseekvationer, som är ett system av andra gradens ordinära dierentialekvationer, vilka generellt är ickelinjära [4].
Ett alternativ till Newtons metod för att ta fram ett systems rörelseek-vationer, att se på framställningen av dynamiken för ett mekaniskt system, är Lagranges ekvationer. Newtons metod baseras på sambandet mellan vek-torstorheter så som momentekvationen som ger oss sambandet mellan vin-kelaccelerationen och kraftmomentet. Istället för att betrakta dessa vektor-storheter använder sig Lagranges metod av systemets energi [8].
Vid en första anblick kan Lagranges formulering av mekaniken se ut som en abstrakt omskrivning av Newtons andra lag, men fördelen i denna ligger i att vilken typ av koordinatsystem vi än väljer att använda så kommer formu-leringen inte att ändras. Newtons vektorstorheter skulle däremot fallerat om det mekaniska systemet istället skulle uttryckas i mer generella koordinater [8]. Då det nns så kallade tvångskrafter, som ofta är okända, förenklar Lag-ranges ekvationer våra beräkningar genom att inte ta hänsyn till dessa. Om systemet skulle bestå av tvångskrafter väljs istället lämpliga generaliserade koordinater och därefter formuleras Lagranges ekvation i termer av de val-da koordinaterna. Metoden ger däremot ingen garanti för att ekvationerna, i förhållande till Newtons, skulle bli enklare att i slutänden beräkna utan fördelen ligger alltså i reduceringen av tvångskrafterna [16].
2.1 Härledning av Lagranges ekvation genom Newtons
andra lag
För att visa att Lagranges ekvationer är oberoende av vilka koordinater som väljs, kommer vi med hjälp av [8] att göra en härledning av Lagranges ekvationer genom Newtons andra lag.
Antalet frihetsgrader för ett system anger hur många oberoende variabler som krävs för att beskriva en kropps position. En partikel som får röra sig fritt i ett tredimenssionellt rum har tre frihetsgrader. För ett partikelsystem som är uppbyggt av N stycken partiklar gäller det att antalet frihetsgrader är 3N. Fortsättningvis kommer antalet frihetsgrader att vara ekvivalent med antalet generaliserade koordinater som krävs för att beskriva ett systems rörelse [16].
Låt ett partikelsystem ha N partiklar, där den n:te partikeln har läget rn= [xn, yn, zn]T. För att lättare kunna förklara förändringar av rörelser för
partikelsystemet vill vi koppla samman dessa kartesiska koordinater med ett generaliserat s-koordinatsystem. Detta gör vi genom att de niera
s1 = x1, s2= y1, s3= z1, ..., sD−2= xN, sD−1= yN, sD = zN, (2.1)
där D = 3N är partikelsystemets frihetsgrader. Sätt s = (s1, s2, s3, ..., sD).
De generaliserade krafterna Fi, för i = 1, ..., D, tillhörande s-koordinatsystemet
kopplas sedan samman med partiklarnas krafter som
F1 = fx1, F2= fy1, F3= fz1, ..., FD = fzN. (2.2)
Även de generaliserade massorna Mi, för i = 1, ..., D, tillhörande s-koordinatsystemet
kopplas samman med massorna för varje partikel, genom
M3K+i= mk+1, (2.3)
där k = 0, 1, ..., N − 1 och i = 1, 2, 3. Med dessa omskrivningar kan Newtons andra lag skrivas om till s-systemet, som
Fi = Mi
d2si
dt2 (2.4)
för i = 1, ..., D.
Den totala kraften fk = [fxk, fyk, fzk]
Tsom verkar på k:te partikeln i ett
partikelsystem, i ett konservativt kraftfält, kan beskrivas genom en partiell derivering av en potentialfunktion, V (r1, ...,rN, t). Dock kan det nnas delar
i fksom inte uppkommer då vi deriverar potentialen, dessa betecknas f (N P ) k . Vi får att fk= − ∂ ∂rk V (r1, ..., rN, t) + f(N P )k . (2.5)
Potentialen V (r1, ...,rN, t)skrivs om till s-systemet genom (2.1) enligt
V (r1, ...,rN, t) = V (s, t) . (2.6)
Genom att använda (2.1), (2.2) och (2.6) fås (2.5) utryckt i s-systemet Fi= −
∂ ∂si
V (s, t) + Fi(N P ) (2.7) för i = 1, ..., D. Från (1.2) kan den kinetiska energin för ett partikelsystem skrivas om till s-systemet som
T = 1 2 D X i=1 Mis˙i2. (2.8)
Utifrån (2.8) och (2.6) de nieras Lagrangefunktionen uttryckt i s-systemet som L (s, ˙s, t) = T ( ˙s) − V (s, t) = 1 2 D X i=1 Mis˙i2− V (s, t) . (2.9)
Genom att partiellt derivera Lagrangefunktionen med avseende på si och ˙si
fås ∂ ∂si L (s, ˙s, t) = − ∂ ∂si V (s, t) (2.10) och ∂ ∂ ˙si L (s, ˙s, t) = ∂ ∂ ˙si 1 2 D X 1 Mis˙i2 = Mis˙i. (2.11)
(2.7) sätts lika med (2.4) och genom att använda (2.10) och (2.11) får vi d dt ∂L(s, ˙s, t) ∂ ˙si −∂L(s, ˙s, t) ∂si = Fi(N P ) (2.12) vilket är Lagranges ekvationer i s-systemet.
Som vi har visat är s-koordinatsystemet en omskrivning av de kartesiska koordinaterna. Vi inför nu de generaliserade koordinaterna, qi. Antag att
varje si kan skrivas som en funktion av q-variabler och tiden t. Sätt q =
(q1, q2, ..., qD). Detta ger
si = si(q, t) (2.13)
för i = 1, ..., D. Vi antar att s(q, t) är inverterbar med avseende på q och ˙q och därför kan q skrivas som en funktion av s, dvs q(s, t).
Vi kan även skriva varje ˙si som
för i = 1, ..., D. Genom att använda ekvationerna (2.13) och (2.14) i (2.9) fås Lagrangefunktionen att bero på q, ˙q och t enligt
L(q, ˙q, t) = L(s(q, t), ˙s(q, ˙q, t), t). (2.15) Genom att använda kedjeregeln på (2.13), med avseende på tiden t, fås att
˙si kan uttryckas som en funktion av q och dess tidsderivator
˙si = dsi(q, t) dt = D X j=1 ∂si(q, t) ∂qj dqj dt + ∂si(q, t) ∂t . (2.16) Låt Q(N P )
j beteckna kraften för ett generaliserat q - koordinatsystem.
Kraften ges av ekvationen
Q(N P )j = D X i=1 Fi(N P )∂si(q, t) ∂qj (2.17)
som har inversen
Fi(N P )= D X j=1 Q(N P )j ∂qj(s, t) ∂si . (2.18)
Lemma 2.1.1. Vi har att
∂ ˙si(q, ˙q, t) ∂ ˙qj = ∂si(q, t) ∂qj och ∂ ˙si(q, ˙q, t) ∂qj = d dt ∂si(q, t) ∂qj . (2.19)
Bevis. Eftersom funktionerna ∂si(q, t)/∂qj och ∂si(q, t)/∂ti ekvation (2.16)
endast beror på q och t innebär det att koe cienten framför ˙qj i (2.16) är
den enda som kvarstår då vi partiellt deriverat med avseende på ˙qj, vilket i
sin tur ger oss bevis på att den första delen av (2.19) håller. För att bevisa den andra delen av (2.19) skrivs VL om som
∂ ˙si(q, ˙q, t) ∂qj = D X k=1 ∂ ∂qj ∂si(q, t) ∂qk dqk dt + ∂ ∂qj ∂si(q, t) ∂t . (2.20)
HL kan, för en godtycklig funktion h(q, t), skrivas om som dh(q, t) dt = D X k=1 ∂h(q, t) ∂qk dqk dt + ∂h(q, t) ∂t . (2.21)
Låt sedan h(q, t) = ∂si(q, t)/∂qj vilket ger d dt ∂si(q, t) ∂qj = D X k=1 ∂ ∂qk ∂si(q, t) ∂qj dqk dt + ∂ ∂t ∂si(q, t) ∂qj (2.22)
det vill säga lika med VL.
2.1.1 Satsen om generellt koordinatbyte
Vi vill med följande sats visa att (2.12) är oberoende av vilka koordinater som väljs.
Sats 2.1.1. Antag att vi gjort ett byte från s koordinatsystemet till q -koordinatsystemet utifrån ekvation (2.13) och de niera Lagrangefunktionen genom (2.15). De niera även den generaliserade kraften för q - koordinatsy-stemet genom ekvation (2.17). Då gäller
d dt ∂L(s, ˙s, t) ∂ ˙si −∂L(s, ˙s, t) ∂si = Fi(N P ) (2.23) om och endast om d dt ∂L(q, ˙q, t) ∂ ˙qk −∂L(q, ˙q, t) ∂qk = Q(N P )k (2.24) för i, j = 1, ..., D.
Bevis. Vi visar först att ekvation (2.23) ger (2.24). Multiplicera både HL och VL i ekvation (2.23) med PD
i=1∂si(q, t)/∂qk, vilket ger D X i=1 ∂si(q, t) ∂qk d dt ∂L(s, ˙s, t) ∂ ˙si − D X i=1 ∂si(q, t) ∂qk ∂L(s, ˙s, t) ∂ ˙si = D X i=1 ∂si(q, t) ∂qk Fi(N P ). (2.25) Låt f och g vara två godtyckliga funktioner. Från produktregeln fås att f (dg/dt) = d(f g)/dt−g(df /dt). Sätt f = ∂si(q, t)/∂qkoch g = ∂L(s, ˙s, t)/∂ ˙si.
Detta ger att den första termen i (2.25) kan skrivas
D X i=1 ∂si(q, t) ∂qk d dt ∂L(s, ˙s, t) ∂ ˙si = D X i=1 d dt ∂si(q, t) ∂qk ∂L(s, ˙s, t) ∂ ˙si − D X i=1 ∂L(s, ˙s, t) ∂ ˙si d dt ∂si(q, t) ∂qk .
Här ersätts både den första och den andra termen i HL av termerna i (2.19) D X i=1 d dt ∂si(q, t) ∂qk ∂L(s, ˙s, t) ∂ ˙si − D X i=1 ∂L(s, ˙s, t) ∂ ˙si d dt ∂si(q, t) ∂qk = D X i=1 d dt ∂ ˙si(q, t) ∂ ˙qk ∂L(s, ˙s, t) ∂ ˙si − D X i=1 ∂L(s, ˙s, t) ∂ ˙si ∂ ˙si(q, ˙q, t) ∂qk . (2.26)
Vidare fås att (2.25) ges av
D X i=1 d dt ∂ ˙si(q, t) ∂ ˙qk ∂L(s, ˙s, t) ∂ ˙si − D X i=1 ∂L(s, ˙s, t) ∂ ˙si ∂ ˙si(q, ˙q, t) ∂qk − D X i=1 ∂si(q, t) ∂qk ∂L(s, ˙s, t) ∂si = D X i=1 ∂si(q, t) ∂qk Fi(N P ) vilket ger d dt D X i=1 ∂L(s, ˙s, t) ∂ ˙si ∂ ˙si(q, t) ∂ ˙qk ! − D X i=1 ∂L(s, ˙s, t) ∂si ∂si(q, t) ∂qk + D X i=1 ∂L(s, ˙s, t) ∂ ˙si ∂ ˙si(q, ˙q, t) ∂qk ! = Q(N P )k . (2.27) Här har HL i (2.27) ersatts av (2.17). Den första parentesen i VL i (2.27) är en utvidgning av kedjeregeln för ∂L(s(q, t), ˙s(q, ˙q, t), t)/∂ ˙qk. Den
and-ra parentesen i (2.27) är en utvidgning av kedjeregeln med avseende på ∂L(s(q, t), ˙s(q, ˙q, t), t)/∂qk. Ovanstående ekvation och (2.15) ger oss
slut-ligen (2.24) d dt ∂L(q, ˙q, t) ∂ ˙qk −∂L(q, ˙q, t) ∂qk = Q(N P )k .
Vi visar nu att (2.24) ger (2.23) genom att göra ovanstående bevis baklänges och vi får d dt ∂L(q, ˙q, t) ∂ ˙qk −∂L(q, ˙q, t) ∂qk = Q(N P )k vilket ger d dt D X j=1 ∂L(s, ˙s, t) ∂ ˙sj ∂ ˙sj(q, t) ∂ ˙qk − D X j=1 ∂L(s, ˙s, t) ∂sj ∂sj(q, t) ∂qk + D X j=1 ∂L(s, ˙s, t) ∂ ˙sj ∂ ˙sj(q, ˙q, t) ∂qk = Q (N P ) k
och D X j=1 d dt ∂ ˙sj(q, t) ∂ ˙qk ∂L(s, ˙s, t) ∂ ˙sj − D X j=1 ∂L(s, ˙s, t) ∂ ˙sj ∂ ˙sj(q, ˙q, t) ∂qk − D X j=1 ∂sj(q, t) ∂qk ∂L(s, ˙s, t) ∂sj = D X j=1 ∂sj(q, t) ∂qk Fj(N P ). (2.28)
Utifrån omskrivningen av ekvation (2.26) gäller
D X j=1 ∂sj(q, t) ∂qk d dt ∂L(s, ˙s, t) ∂ ˙sj − D X j=1 ∂sj(q, t) ∂qk ∂L(s, ˙s, t) ∂ ˙sj = D X j=1 ∂sj(q, t) ∂qk Fj(N P ). (2.29)
Slutligen multipliceras (2.25) med inversen till ∂sj(q, t)/∂qk, vilket är ∂qk(s, t)/∂sj
där k = 1, ..., D, och vi har d dt ∂L(s, ˙s, t) ∂ ˙sj −∂L(s, ˙s, t) ∂sj = Fj(N P ).
2.2 Hamiltons princip
Vi har, med hjälp av Newtons andra lag, visat att Lagranges ekvationer är oberoende av vilka generaliserade koordinater som väljs. En mer gene-rell metod för att nna rörelseekvationer är att använda sig av Hamiltons princip, vilken vi utifrån [4] kommer att härleda. Hamiltons princip är en grundläggande princip inom mekaniken som beskriver ett systems rörelse från punkt A till punkt B, mellan tidpunkterna t0 och t1. Genom den så
kallade verkningsintegralen
S = Z t1
t0
L(q, ˙q, t) dt fås Hamiltons variationsprincip enligt
δS = δ Z t1
t0
L(q, ˙q, t) dt = 0 (2.30) där δS är variationen av S. Det q som uppfyller (2.30) ger den rörelse mellan Aoch B där S har ett extremum, vilket också benämns den verkliga vägen.
Antag att variationen av en generaliserad koordinat qi kan skrivas som
en funktion av tiden t på intervallet [t0, t1]och där variationen av denna kan
skrivas enligt
δqi(t) = ηi(t) (2.31)
där ηi(t)är en godtycklig kontinuerlig funktion som uppfyller ηi(t0) = ηi(t1) =
0och är godtyckligt liten. Antag att variationen av S, för något q, kan skri-vas som δS = Z t1 t0 ∂L ∂qη + ∂L ∂ ˙qη˙ dt. (2.32)
Genom att använda (2.32) samt omskrivningen i (2.31) och därmed beräkna variationen av S fås att δS = δ Z t1 t0 L(q, ˙q, t) dt = Z t1 t0 X i ∂L ∂qi δqi+ ∂L ∂ ˙qi δ ˙qi dt = 0 (2.33) där vi enligt (2.31) får att δ ˙qi = d dtηi(t) = d dtηi(t) = d dtδqi. (2.34)
Vidare sätts (2.34) in i den andra termen i (2.33) och vi får genom partiell integration av den termen att
Z t1 t0 X i ∂L ∂ ˙qi δ ˙qidt = Z t1 t0 X i ∂L ∂ ˙qi d dt(δqi) dt = " X i ∂L ∂ ˙qi δqi #t1 t0 − Z t1 t0 X i d dt ∂L ∂ ˙qi δqidt. (2.35)
Eftersom att δqi(t0) = δqi(t1) = 0kommer uttrycket inom klammern i (2.35)
att bli noll. Det gäller då att δS = δ Z t1 t0 L(q, ˙q, t) dt = Z t1 t0 X i ∂L ∂qi − d dt ∂L ∂ ˙qi δqidt = 0. (2.36)
De generaliserade koordinaterna är oberoende och således även variatio-nerna av dem. På så sätt kan till exempel δqk 6= 0 och δqi = 0, i 6= k på
[t0, t1], vilket medför i (2.36) att
Z t1 t0 ∂L ∂qk − d dt ∂L ∂ ˙qk δqkdt = 0. (2.37)
Då vi vet att δqk 6= 0 kan väljas godtyckligt innebär det att klammern i
(2.37) måste vara noll på hela [t0, t1]. På samma sätt fås att koe cienterna
framför δq2, δq3, ..., δqs måste vara noll, dvs att
∂L ∂qi − d dt ∂L ∂ ˙qi = 0 (2.38)
för i = 1, ..., s.
Lösningen till Lagranges ekvation ger oss s ordinära dierentialekvatio-ner, vilket beskriver systemets rörelse längs med den verkliga vägen. Den-na rörelse motsvarar kraften som uppkommer från den potentiella energin i Lagrangefunktionen. Det som ovanstående härledning visar är alltså att δS = 0, gäller endast då variationerna av en vald väg motsvaras av den verkliga vägen. Detta ekvationssystem kallas för Lagranges ekvationer för ett konservativt system.
Vi vill belysa att ovanstående bevis håller i vårt fall, då endast kontinu-erliga funktioner används. I annat fall kan det nnas diskreta funktioner då beviset inte längre håller.
Här ses även att termerna i VL i (2.38) har bytt plats i förhållande till tidigare härledning. Detta beror på olika traditioner att använda sig av Lag-ranges ekvationer. Inom mekaniken används formuleringen enligt härledning-en från Newtons andra lag medan dhärledning-en inom variationskalkylhärledning-en formuleras enligt härledningen från Hamiltons princip.
Det vi även kan notera är att utifrån vår härledning från Hamiltons prin-cip är högerledet noll på grund av att det är ett system utan geometriska tvångsvillkor; villkor som innebär att en partikel (eller kropp) endast får röra sig på ett visst sätt, på en viss yta eller liknande. Eventuella tvångsvillkor delas in i holonoma och icke-holonoma tvång, där de holonoma tvången in-nebär ett samband mellan systemets koordinater och tiden [4]. Om systemet skulle påverkas av så kallade icke-holonoma tvångsvillkor, det vill säga att det nns ett samband mellan koordinaterna och tiden men även med syste-mets hastigheter, används istället härledningen enligt [10], där högerledet blir lika med någon generaliserad kraft Qj.
2.3 Ett exempel på Lagranges metod
Antag att vi har ett tvådimensionellt koordinatsystem. Låt en stång, med längden l, fästas i Origo så att den hänger i negativt y - led och låt ϕ vara vinkeln mellan stången och den positiva x - axeln. Masscentrum, G, antas sitta i mitten av stången (se gur 2.1).
Vi vill använda Lagranges funktion, L = T − V , för att sedan härleda systemets rörelseekvationer. Vi behöver därför den kinetiska samt den poten-tiella energin, T respektive V . Enligt de nitionen för kinetisk energi gäller att T = Iω2/2, där ω = ˙ϕ är systemets vinkelhastighet. För den potentiella
energin gäller enligt de nition att V = mgyG, där yG är stångens
masscent-rums y-koordinat. Enligt de nitionen för tröghetsmoment I för en stång som är upphängd i sin ena ände fås I = ml2/3, se [4]. Vi får alltså
L = T −V = 1 2Iω 2−mgy G= 1 2 ml2 3 ϕ˙ 2−(−mgl 2sin ϕ) = ml2 6 ϕ˙ 2+mgl 2sin ϕ.
Figur 2.1: Pendlande stång. ϕ = vinkel mellan stången och den positiva x -axeln, l = stångens längd
För att ta reda på rörelseekvationerna för systemet sätts ovanstående funk-tion in i Lagranges ekvafunk-tion. Vi antar Q(N P )
j = 0 och får d dt ∂L(ϕ, ˙ϕ, t) ∂ ˙ϕ −∂L(ϕ, ˙ϕ, t) ∂ϕ = Q (N P ) j vilket ger d dt ml2 3 ϕ˙ − mgl cos ϕ = 0 och vi får ml2 3 ϕ − mg¨ l 2cos ϕ = 0.
För en entydig lösning till denna dierentialekvation krävs även två begyn-nelsevillkor. Om vi nu antar att begynnelsevillkoren för stången är ˙ϕ(0) = 0 samt ϕ(0) = π/4, får vi en lösning på rörelseekvationerna. Vi löser ut ¨ϕoch får
¨
ϕ = 3g cos ϕ 2l .
Sedan inför vi två nya variabler, u1 och u2, och sätter
u1 = ϕ, u2= ˙ϕ vilket ger ˙ u1 = u2 = ˙ϕ, ˙u2= ¨u1= ¨ϕ = 3g cos ϕ 2l = 3g cos u1 2l .
Detta ger oss istället för en andra ordningens dierentialekvation, två stycken första ordningens dierentialekvationer som går att lösa med lösarenode45
Figur 2.2: Position för stångens masscentrum vid olika tillfällen under rörel-sen med begynnelsevillkoren ϕ(0) = π/4 och ˙ϕ(0) = 0.
1 [TOUT,YOUT] = ode45('fpendelum',[0; 1],[pi/4; 0]),
där inte alla kommandorader är utskrivna. Det som då behövs är en Mat-labfunktionfpendelum(t,yin)som är just dessa två dierentialekvationer,
tidsintervallet och de begynnelsevillkor vi beskrivit ovan. Vi har antagit att l = 0.01och g = 9.81 och vi plottar sedan ut positionen för stångens mass-centrum under rörelsen. Fullständig Matlabkod kan ni se i bilaga A.
Om begynnelsevillkoren ändras så att till exempel ϕ(0) = 0, resuluterar det i att rörelsen för stången blir längre och att stången därmed inte hinner pendla lika många gånger under tidsintervallet.
Figur 2.3: Vinklarna ϕ och ˙ϕ under rörelsen, med begynnelsevillkoren ϕ(0) = π/4och ˙ϕ(0) = 0.
Figur 2.4: Position för stångens masscentrum vid olika tillfällen under rörel-sen med begynnelsevillkoren ϕ(0) = 0 och ˙ϕ(0) = 0.
Figur 2.5: Vinklarna ϕ och ˙ϕ under rörelsen, med begynnelsevillkoren ϕ(0) = 0och ˙ϕ(0) = 0.
Kapitel 3
Släggkastets rotationsfas
utifrån Lagranges metod
För att ta fram en modell som bäst beskriver släggkastets rotationsfas, men begränsa ramen för arbetet, har vi valt att utgå ifrån stela kroppar. Män-niskokroppen som håller i släggan får i denna uppsats antas utgöras av en cylinder medan slägghuvudet antas vara en punktmassa fäst i en masslös tråd.
Anledningen till att vi inte använder oss av människokroppen som den ser ut i verkligheten är att vi bland annat då måste ta hänsyn till muskel-och senkrafter samt yttre krafter som kan komma att forma eller påverka kroppen på annat sätt. Det skulle samtidigt innebära att tröghetsmomentet för kroppen blir alltför kompicerat att beräkna. Att slägghuvudet antas vara fäst i en masslös tråd beror på att beräkningen då kan göras på tröghetsmo-mentet som uppstår kring punktmassan samt kring cylindern. Vi kan alltså bortse från tröghetsmomentet kring vajern som släggan vanligtvis är fäst vid, vilket ger oss enklare beräkningar. Vi kommer i beräkningarna inte att ta hänsyn till yttre krafter så som friktion och luftmotstånd.
Vi antar att de så kallade försvängarna i släggkastet redan gjorts så att punktmassan roterar med en viss hastighet vilket ger begynnelsevärden i modellen.
3.1 Beräkning av rörelseekvationer för vår modell
Antag att vi har en homogen cylinder C med masscentrum GC, vars ena ändeär fäst i Origo. Cylindern antas ha radien r, höjden l1, massan m1 och
voly-men V . Vi väljer ett koordinatsystem OXY Z samt Ox0y0z0 med gemensamt
origo och Ox00y00z00, så att Ox0y0z0är xt i cylindern och Ox00y00z00 xt i
punkt-massan, (se gur 3.1). Sätt c1 = cos(θ1), c2 = cos(θ2), c3 = cos(θ3), c4 =
cos(θ4), s1 = sin(θ1), s2 = sin(θ2), s3= sin(θ3) och s4= sin(θ4) där θ1, θ2, θ3
Figur 3.1: Modell av släggkastaren. För cylindern gäller: r = radien, l1 =
höjden, m1 = massan, V = volymen och GC = masscentrum. För
punkt-massan gäller: m2=massan, GP =masscentrum och l2 =längden av tråden
ut till punktmassan.
Figur 3.2: Modell av släggkastaren utifrån Ox0y0z0. För cylindern gäller: l 1 =
höjden och m1 = massan. För punktmassan gäller: m2 = massan och l2 =
längden av tråden ut till punktmassan. Vinkeln mellan z0 - axeln och
punkt-massan de nieras θ4. x0 = l2sin θ4 och z0 = l1+ l2cos θ4 är det vinkelräta
θ4 införs för att beskriva punktmassans position i förhållande till Ox0y0z0.
Utifrån Eulervinklarna fås vinkelhastighetsvektorn ω0 enligt (1.8)
ω0 = ωx ωy ωz = ˙ θ1s2s3+ ˙θ2c3 ˙ θ1c3s2− ˙θ2s3 ˙ θ1c2+ ˙θ3 .
På toppen av cylindern fästs en tråd, med längden l2, i ena änden och där
det i den andra änden av tråden fästs en punktmassa, P , med massan m2och
masscentrum GP. Punktmassans ursprungsläge antas beskrivas genom
posi-tionsvektorn r0 = [x0, 0, z0]T. Fortsättningsvis, då vi pratar om vår kropp,
menar vi cylindern och punktmassan (inkl. den masslösa tråden).
Låt nu kroppen rotera likt det tidigare beskrivna släggkastet. Vinkeln θ1
beskriver rotationen kring Z - axeln. Den vinkel som beskriver rotationen i fötterna/benen på kastaren antas vara θ3. Vinklarna θ2 och θ4 beskriver
lutningen på kastaren i förhållande till Z - axeln respektive vinkeln mellan z0 - axeln och slägghuvudet. Släggkastarens rörelse mot kastriktningen, som vi antar är i negativt x - led, kommer att beskrivas med hjälp av en transla-torisk rörelse där v = [−vX, 0, 0]T antas vara konstant. Beräkningar av den
kinetiska energin görs utifrån det xa koordinatsystemet Ox0y0z0.
3.1.1 Cylinderns kinetiska energi
Utifrån (1.1) fås rotationsenergin som TrotC = ω
0TI0
Cω02/2 där
tröghetsmo-menten kring respektive axel ges av Ix0x0 = Z y2+ z2dm = Z Z Z C y2+ z2 ρ dx dy dz = m1 πr2l 1 Z Z Z C y2+ z2dx dy dz = m1 πr2l 1 πl1 r4 4 + πr 2l31 3 = m1 12 3r 2+ 4l2 1 .
Tack vare symmetri kring x0 - respektive y0- axeln är I
y0y0 = Ix0x0. Vidare fås att Iz0z0 = Z x2+ y2dm = Z Z Z C x2+ y2 ρ dx dy dz = m1 πr2l 1 Z Z Z C x2+ y2dx dy dz = m1 πr2l 1 πr4l1 2 = m1r 2 2 .
Efter ytterligare beräkningar fås att tröghetsprodukterna Ix0y0 = Ix0z0 = Iy0x0 = Iy0z0 = Iz0x0 = Iz0y0 = 0 vilket ger oss tröghetstensorn
I0C = Ix0x0 0 0 0 Iy0y0 0 0 0 Iz0z0 .
Detta ger TrotC = 1 2Ix0x0 ωx20 + ωy20 +12Iz0z0(ωz0)2 = 12Ix0x0 ˙θ2 2+ ˙θ21sin2(θ2) +12Iz0z0 ˙θ3+ ˙θ1cos(θ2) 2 där Ix0x0, Iz0z0 är oberoende av θi.
Utifrån (1.1) samt antagandet att v är konstant ges translationsenergin uti-från TtraC = 1 2(m1+ m2) v Tv = 1 2(m1+ m2) v 2 X
så att cylinders totala kinetiska energi är TC = TtraC+ TrotC = 12(m1+ m2) vX2 +12Ix0x0 ˙θ22+ ˙θ21sin2(θ2) +12Iz0z0 ˙θ3+ ˙θ1cos(θ2) 2 . (3.1)
3.1.2 Punktmassans kinetiska energi
Punktmassans hastighetsvektor, v0 =v
x0, vy0, vz0T, i Ox0y0z0 - systemet ges dels av rotationen kring toppen av cylindern, dels av kroppens rotation, dvs w0 uttryckt i Ox0y0z0 - systemet. Vi får då att den totala hastigheten ges av (1.3) så att v0 = v0+ w0× r0 där v0 = ˙ x0 0 ˙ z0 , r0 = x0 0 z0
och där x0= l2sin(θ4) och z0= l1+ l2cos(θ4)samt att ˙x0= l2cos(θ4) ˙θ4 och
˙
z0 = −l2sin(θ4) ˙θ4. Detta ger
vx0 = ˙x0+ wy0z0 = l2cos(θ4) ˙θ4+ ˙θ1c3s2− ˙θ2s3 (l1+ l2cos(θ4)) , vy0 = wz0x0− wx0z0 = ˙θ1c2+ ˙θ3 l2sin(θ4) − ˙θ1s2s3+ ˙θ2c3 (l1+ l2cos(θ4)) , vz0 = ˙z0− ωy0x0 = −l2sin(θ4) ˙θ4− ˙θ1c3s2− ˙θ2s3 l2sin(θ4).
Vi får utifrån (1.2) att den kinetiska energin för punktmassan ges av TP = 12m2v0Tv0 = 12m2 v2 x0 + vy20+ vz20 = 12m2 x˙0+ ωy0z02+1 2m2(ωz0x0− ωx0z0) 2 +12m2 z˙0− ωy0x02 vilket kan förenklas till
TP = 12m2 x˙20+ ˙z02 + m2 x˙0wy0z0− ˙z0wy0x0 +12m2 w2y0z20+ (wz0x0− wx0z0)2+ w2 y0x20 . (3.2)
3.1.3 Cylindern och punktmassans potentiella energi
Utifrån (1.4) fås att den potentiella energin hos cylindern VC och
punktmas-san VP är VC = mgZG= m1g l1 2 cos θ2 (3.3) respektive VP = mgZG= m2g(l1cos θ2+ l2cos θ4). (3.4) 3.1.4 Lagrangefunktionen för kroppen
För att ta fram den fullständiga Lagrangefunktionen L = T − V , används T = TC + TP från (3.1) samt (3.2) och V = VC+ VP från (3.3) samt (3.4),
vilket ger L = T − V = TC+ TP − (VC + VP) = 12Ix0x0 ωx20 + ωy20 +12Iz0z0(ωz0)2 +12m2l22θ˙42+ m2 x˙0wy0z0− ˙z0wy0x0 +12m2 w2y0z20+ (wz0x0− wx0z0)2+ w2 y0x20 = 12θ˙12 s22 Ix0x0 + m2l21 + c22Iz0z0+ c24s22s23m2l22+ c23s22m2l22 +c22s24m2l22+ 2c4s22m2l1l2− 2c2s2s3s4m2l1l2− 2c2c4s2s3s4m2l22 +12θ˙22 Ix0x0+ m2l2 1+ c23c24m2l22+ s23m2l22+ 2c4m2l1l2 +12θ˙23 Iz0z0+ s2 4m2l22 +12θ˙24 m2l22 + ˙θ1θ˙2 c3c24s2s3m2l22− c2c3s4m2l1l2− c3s2s3m2l22− c2c3c4s4m2l22 + ˙θ1θ˙3 c2Iz0z0+ c2s2 4m2l22− s2s3s4m2l1l2− c4s2s3s4m2l22 + ˙θ1θ˙4 c3c4s2m2l1l2+ c3s2m2l22 − ˙θ2θ˙3 c3s4m2l1l2+ c3c4s4m2l22 − ˙θ2θ˙4 c4s3m2l1l2+ s3m2l22 −m1gl21c2+ m2g (l1c2+ l2c4) . 3.1.5 Rörelseekvationer
Vi antar att vinklarna θ1, θ2, θ3 och θ4 är generaliserade koordinater där de
tagit fram Lagrangefunktionen för kroppen vill vi ta fram Lagranges ekvatio-ner, det vill säga kroppens rörelseekvationer. För att göra detta bestäms de partiella derivatorna av L, för varje generaliserad koordinat q1, q2, q3 och q4,
vilka vi behöver för att formulera rörelseekvationen för respektie koordinat. Från (2.24) fås rörelseekvationen för respektive generaliserad koordinat.
Genom att först derivera L med avseende på ˙qi och därefter derivera det
vi får ut, med avseende på tiden t fås den första termen i VL i Lagranges ekvation. Den andra termen i VL fås genom att vi deriverar L med avseende på qi. Den generaliserade kraften Q(N P )i sätts som givna moment Mi(t)och
vi får d dt ∂L ∂ ˙qi − ∂L ∂qi = Mi(t)
där i = 1, ..., 4. För resultat av dessa rörelseekvationer se bilaga C.
3.2 Lagranges ekvationer på matrisform
Vår uppsats har genomgående handlat om att ta fram Lagranges ekvationer, från ett systems kinetiska -och potentiella energi, för att få dess rörelseekva-tioner. Det vi efter beräkningar har insett är att arbetet med denna metod och dessa rörelseekvationer blir alltför omfattande för att kunna göra en si-mulering i Matlab. Till följd av detta har vi istället valt att använda oss av en annan metod där Lagranges ekvationer skrivs på matrisform och fortsätt-ningsvis antar vi att θ4 är konstant, dvs ˙θ4= 0.
Vid implementering av Lagranges ekvationer i Matlab är det enklast att de niera ekvationerna med matriser och vektorer. Detta kan göras på ett kompakt och systematiskt sätt, se [1] för en mer detaljerad beskrivning. Här ska vi endast ge de ekvationer som fås.
Antag att a är en vektor
a = a1 a2 a3 och vi inför beteckningen
[a] = 0 −a3 a2 a3 0 −a1 −a2 a1 0 .
Med ω betecknar vi den i det kropps xa systemet givna rotationsvektor som fås då det kropps xa systemet roterar med Eulervinkelhastigheterna och vi har
där W = s2s3 c3 0 s2c3 −s3 0 c2 0 1 , ˙θ = ˙ θ1 ˙ θ2 ˙ θ3
och c1 = cos(θ1), s1= sin(θ1), c2 = cos(θ2), . . .. Antag att tröghetsmomentet
för cylindern är IC och för punktmassan IP taget i det kropps xa
koordi-natsystemet och sätt IT ot = IC+ IP. Då fås följande Lagrangeekvationer
WTIT otW ¨θ + (WTIT otW + W˙ T[W ˙θ]IT otW ) ˙θ − ∂V (θ) ∂θ = Q där ˙ W =X j ∂W ∂θj ˙ θj och V (θ) är potentialen.
I Matlab löses detta, andra ordningens system, som ett system av första ordningens ekvationer genom att införa en ny variabel σ = ˙θ vilket ger
WTIT otW ˙σ + (WTIT otW + W˙ T[W σ]IT otW )σ − ∂V (θ) ∂θ = Q (3.5) och ˙ θ = σ.
3.3 Implementering i Matlab av
matrisformulering-en
I Matlab nns ett ertal lösare till system av dierentialekvationer som t.ex.
ode45som kan hantera systemet (3.5). I Matlabs notation skrivs system som
1 M(t,y)*y' = f(t,y)
däry,M(t,y)och f(t,y)motsvaras av
θ σ , WTItotW 0 0 I , " −(WTI
totW + W˙ T[W σ]ItotW )σ +
∂V ∂θ σ
# . Förutom en Matlabfunktion förf(t,y)krävs en Matlabfunktion förM(t,y)
enligt följande (från help ode45)
1 ode45 can solve problems M(t,y)*y' = f(t,y) with mass ...
matrix M that is nonsingular. Use ODESET to set the ...
'Mass' property to a function handle MASS if MASS(T,Y) ...
Vi behöver alltså använda odeset vilket går till på följande sätt, där inte
alla kommandorader är med
1 options=odeset('mass',@name_of_M_function)
2 [t,y]=ode45(@name_of_f_function,tspan,y0,options)
För fullständig Matlabkod se bilaga B.
3.4 Avancerad simulering
Den modell vi använt är en grov förenkling av den verkliga rörelsen. För att få en mer realistisk modell blir vi tvungna att modellera kroppen mer detal-jerat med hänsyn till leder och muskler. Ett sådant simuleringsverktyg är [3] som har använts i era sammanhang som t.ex. skidstakning med intressanta resultat, se [11],[14].
Ett alterantivt simuleringsverktyg till Matlab är ODE [15]. I ODE ba-seras simuleringen på de s.k. Euler-Newtons ekvationer [17] vilket kort sagt integrerar rörelsen givet krafter och moment (krafter och moment ger via integration translationsacceleration och vinkelaccelreation som i sin tur ger translationshastighet och vinkelhastighet, vilket slutligen ger rörelsen).
Via länken http://youtu.be/HKzoZWa7aXM ses en enkel simulering i ODE av ett släggkast.
3.5 Resultat
För att kunna plotta ut slägghuvudets rörelse i Matlab behöver vi anta värden på de komponenter som används. Vi antar att släggkastaren är en man och att längden från fötterna upp till axlarna på han är 1.70 m, vilket ger längden på cylindern l1 = 1.70. Vi har även antagit att vikten på en
manlig släggkastare är cirka 100 kg och att radien på dennes kropp är 0.20 m. Dessa värden motsvarar i vår modell m1 = 100 respektive r = 0.20 hos
cylindern. Längden på en släggas vajer är 1.20 m och vi har antagit att släggkastarens armar är 0.60 m. Dessa två värdena adderas, vilket ger att l2 = 1.80 m. Vikten på slägghuvudet är 7.26 kg, vilket ger värdet på m2 =
7.26. Vi har även antagit att gravitationskraften g = 9.81 och som vi tidigare angivit antas ˙θ4 = 0 då vi antar att θ4 = π2 under hela rörelsen. I modellen
används en konstant translatorisk rörelse i negativt x-led. Tiden det tar för en släggkastare att rotera 4.5 varv och utföra kastet uppskattas till ungefär 3 sekunder, vilket betyder att t ∈ [0, 3] i våra resultat. Begynnelsevärdena för vinklarna antar vi är θ1(0) = 0, θ2(0) = π6, θ3(0) = 0, ˙θ1(0) = 0, ˙θ2(0) =
0, ˙θ3(0) = 24π.
För att få plotten att likna slägghuvudets rörelse under rotationsfasen i ett släggkast, använder vi oss av momenten, Mi(t) där i = 1, 2, 3, och
Figur 3.3: Vinklarnas värde under rotationsfasen i ett släggkast.
begynnelsevärdet för ˙θ3. Dessa moment och begynnelsevärdet för ˙θ3påverkar
i sin tur vinklarna, θ1, θ2och θ3på olika sätt. Genom att sätta rätt värden på
momenten (de värden på momenten som vi använder oss av går att se i bilaga B) vill vi få θ3 att växa hela tiden, vilket betyder att släggkastaren roterar
kring sin egen axel. För att få slägghuvudet att rotera 4.5 varv behöver θ3
bli 4.5 · 2π ≈ 28, 27 under de 3 sekunder som rotationerna ska ske.
Samtidigt som θ3 ska växa får θ1 och θ2 inte bli för stora. Vi vill hålla
θ1 litet och θ2 nära π/6, då θ1 är släggkastarens rotation kring Z - axeln
och vi vill endast att släggkastarens rotationer ska ske i fötterna och på så sätt kring sin egen axel. Om släggkastaren samtidigt skulle rotera kring Z -axeln genererar det en felaktig rörelse. Att hålla θ2 nära π/6 är viktigt då θ2
är vinkeln mellan släggkastaren och Z - axeln. Vi vill ha en liten lutning på släggkastaren, men om θ2 blir för stor betyder det att han faller mot marken.
Därför har vi anpassat momenten för att få de värdena vi vill ha på θ1, θ2
och θ3 och få den rörelsen på slägghuvudet som vi vill.
I plotten börjar θ3 växa direkt men ändå börjar inte slägghuvudet rotera
direkt i plotten. Vi vet inte varför detta sker men bortser vi från den korta tid som vi bara har en translatorisk rörelse, ser resterande del av plotten ut som rotationsfasen i ett släggkast.
I gur 3.3 presenteras vinklarnas värden under de 3 sekunder som ro-tationen sker. Som önskat växer θ3 under hela intervallet upp mot värdet
28.27. θ1 och θ2 varierar istället lite, men håller sig runt de värden vi vill
ha dem. Figur 3.4 - 3.7 visar slägghuvudets position under de 4.5 rotationer som sker, sett ifrån olika håll. Slägghuvudet rör sig likt rörelsen i ett slägg-kast, från en låg position när släggkastaren har ryggen mot kastriktningen
Figur 3.4: Positionen för slägghuvudet vid olika tillfällen under rotationsfasen i ett släggkast sett snett bakifrån.
Figur 3.5: Positionen för slägghuvudet vid olika tillfällen under rotationsfasen i ett släggkast sett från sidan.
Figur 3.6: Positionen för slägghuvudet vid olika tillfällen under rotationsfasen i ett släggkast sett bakifrån.
Figur 3.7: Positionen för slägghuvudet vid olika tillfällen under rotationsfasen i ett släggkast sett uppifrån.
(negativ x-led) till en högre position när släggkastarens kropp roterar mot kastriktningen. Efter 4 varv roteras släggkastarens kropp det sista halva var-vet mot kastriktningen och släggan rör sig mot en högre position innan den släpps iväg.
Kapitel 4
Diskussion
Bakgrunden till vårt val av ämne går tillbaka till det genuina idrottsintres-se vi har, som vi om möjligt ämnat att koppla samman med matematiken. För oss resulterade det i att vi såg möjligheten att använda släggkastet som något att utgå ifrån. Kontentan av det blev att vi tillägnat oss kunskaper i mekanik för att ta fram stelkroppsmodellen och nna en användbar me-tod för att matematiskt beskriva släggkastets rotationsfas. Ytterligare ett led i processen har varit att tillägna oss kunskaper i matematisk simulering i Matlab för att på så sätt få fram en bild av resultatet. Arbetet har varit omfattande och tidskrävande dels på grund av den ovan beskrivna proces-sen, dels för att komplexiteten i de rörelseekvationer som vi tagit fram med Lagranges metod blev för avancerade att lösa i Matlab. Av den anledningen blev vi tvungna att tänka om och försöka tillämpa en annan metod för att simuleringen skulle bli möjlig. Det vi istället har gjort är att använda oss av Lagranges ekvationer på matrisform och därmed kan sägas att de resultat vi åstadkommit är framtagna med hjälp av implementering av dessa matriser i Matlab.
Inspiration till uppsatsen är som sagt hämtad ifrån vårt intresse för id-rott, men även hämtad från en studie, se [12]. I denna studie används en modell, vilken utgörs av 15 sammanlänkade segment. Detta anser vi går att jämföra med vår modell, men där vi gör modellen mer simpel med endast två segment. En skillnad mellan vår uppsats och studien är att vi använder oss av Lagranges ekvationer medan studien använder sig av Newtons och Eulers ekvationer. Studiens resultat och diskussion utgår ifrån hur rörelsen på mas-scentrum för släggan respektive släggkastaren, påverkar kroppens drivande moment i olika positioner av rörelsen.
Den utvecklingspotential vi ser med vår uppsats utifrån den stelkropps-modell vi tagit fram, är att denna skulle kunna bli mer lik en människokropp, som i [12]. Tanken med att få modellen att gestalta en människokropp i så stor utsträckning som möjligt, är att i slutänden kunna genomföra en analys över tekniken i ett släggkast med hjälp av matematiska beräkningar. I
lik-het med [12] som presenterar till exempel att det är fötternas rotation som påverkar ingången i varje ny rotation, snarare än de drivande momenten i släggkastarens överkropp. Det skulle vara intressant att utifrån detta även ta hänsyn till de yttre krafter som påverkar kroppen, så som luftmotstånd och friktion, som vi i uppsatsen har förbisett. Ytterligare en utvecklingspotential vi ser är att matematiskt istället använda sig av kvaternioner. De är enkla att använda och inga singulariteter nns, därför är de också ett populärt alternativ för att beskriva stela kroppars rörelser [9].
Litteraturförteckning
[1] Advanced classical mechanics/rigid bodies. [Hämtad:2014-01-28]
http://en.wikiversity.org/wiki/Advanced_Classical_Mechanics/Rigid_Bodies
[2] Analytisk mekanik. [Hämtad:2013-12-06]
http://www.ne.se.db.ub.oru.se/analytisk-mekanik.
[3] Anybody Technology. [Hämtad:2014-01-28]
http://www.anybodytech.com/
[4] Apazidiz, N (2012): Mekanik II - Partikelsystem, stel kropp och ana-lytisk mekanik. Lund: Studentlitteratur AB
[5] Borgström, A (1998): Utveckling av tekniken i spjutkastning. Svensk idrottsforskning. Nummer 1, sid. 25-30. [Hämtad: 2012-12-10] http://www.gih.se/Documents/CIF/Tidningen/1998/1-1998/SVIF19981hela.pdf
[6] Bull Andersen, T och Kristensen, L B (2007): Biomekanik och rörel-selära - analys av människans rörelser. Stockholm: Liber AB
[7] Böiers, L-C och Persson, A (2005): Analys i era variabler. Lund: Studentlitteratur AB
[8] Davis Johns, O (2005): Analytical mechanics for relativity and quan-tum mechanics. New York: Oxford University press
[9] Diebel, J (2006): Representing attitude: Euler angles, unit quater-nions and rotation vectors. Stanford, California: Stanford Univer-sity. [Hämtad: 2014-01-15] ftp://sbai2009.ene.unb.br/Projects/GPS-IMU/George/arquivos/Bibliogra a/79.pdf
[10] Goldstein, H, Poole, C och Safko, J (2002): Classical mechanics. Up-per Saddle River: Pearson education.
[11] Holmberg, L J och Lund, A M (2008): A musculoskeletal full body simulation of crosscountry skiing. Proc. IMechE, Part P: J. Sports Engineering and Technology. Volym 222, upplaga 1, sid. 11-22
[12] T. K. Karalis, X (1991): Control torque components of center-of-mass motions in hammer throwing. Archive of Applied Mechanics. Volym 61, upplaga 5, sid. 344-360
[13] Li, Z, Murray, R M och Sastry, S S (1994): Mathematical introduction to robotic manipulation. Pasadena, California: California institute of technology. [Hämtad: 2014-06-09] http://www.cds.caltech.edu/ mur-ray/mlswiki
[14] Lund, M, Ståhl, F och Gulliksson, M (2008): Regularity aspects in inverse musculoskeletal biomechanics. Numerical analysis and applied mathematics, (AIP conference proceedings). sid. 368-371
[15] Open dynamics engine. [Hämtad:2014-05-07] http://ode.org/
[16] Taylor, J R (2005): Classical mechanics. Sausalito, California: Uni-versity science books
[17] Winter, D A (2009): Biomechanics and motor control of human mo-vement. Hoboken, New Jresey: John Wiley and sons
Bilaga A
Matlabkod för ett exempel på
Lagrange metod
För Matlabfunktionenfpendelumär kommandoraderna
1 function yout = fpendelum(t,yin) 2 g = 9.81; 3 L = 0.01; 4 const = 3*g/(2*L); 5 yout = zeros(2,1); 6 yout(1) = yin(2); 7 yout(2) = const*cos(yin(1)); 8 end
Sedan för att plotta ut rörelsen används kommandoraderna
1 clear 2 format short e 3 clf 4 5 L = 0.01; 6 7 [TOUT,YOUT] = ode45('fpendelum',[0; 1],[0; 0]); 8 9 plot(TOUT,YOUT(:,1),TOUT,YOUT(:,2)) 10 11 y = −L*sin(YOUT(:,1)); 12 x = L*cos(YOUT(:,1)); 13 14 axis([−0.02 0.02 −0.02 0.02]) 15 hold on 16 n = length(TOUT); 17 for k=1:n 18 pause(0.1) 19 plot(x(k),y(k),'*') 20 end
21 hold off
Där man använder plot(TOUT, YOUT(:,1), TOUT, YOUT(:,2)) för att
plotta ut vinklarna och använder det som står under det för att plotta ut rörelsen.
Bilaga B
Matlabkod för släggkastets
rotationsfas
För att plotta ut slägghuvudets rörelse i släggkastets rotationsfas behövs Matlabfunktionerna M(t,y),f(t,y). Kommandoraderna förM(t,y) är
1 function M = lhs_hammer_anton(t,y,I_tot) 2 3 for i=5:6 4 c(i−3)=cos(y(i)); 5 s(i−3)=sin(y(i)); 6 end 7 W=[s(2)*s(3) c(3) 0; c(3)*s(2) −s(3) 0; c(2) 0 1]; 8
9 M=[W'*I_tot*W zeros(3,3); zeros(3,3) eye(3,3)]; 10 11 end Kommandoraderna förf(t,y)är 1 function f = rhs_hammer_anton(t,y,I_tot) 2 3 g=9.81; m_1=100; m_2=7.26; L1=1.70; 4 dV=−sin(y(5))*(m_1*g*L1/2+m_2*g*L1); 5 M1 = 720; 6 M2 = 650; 7 M3 = 750; 8 9 10 for i=5:6 11 c(i−3)=cos(y(i)); 12 s(i−3)=sin(y(i)); 13 end 14 W=[s(2)*s(3) c(3) 0; c(3)*s(2) −s(3) 0; c(2) 0 1]; 15 16 dW3=[s(2)*c(3) −s(3) 0; −s(2)*s(3) −c(3) 0; 0 0 0];
17 dW2=[c(2)*s(3) 0 0; c(2)*c(3) 0 0; −s(2) 0 0]; 18 19 w=W*y(1:3); 20 G=[0 −w(3) w(2); w(3) 0 −w(1); −w(2) w(1) 0]; 21 fvec = −(W'*I_tot*(dW3*y(3)+dW2*y(2)) + ... W'*G*I_tot*W)*y(1:3)+dV; 22 fvec1 = fvec + [M1,M2,M3]'; 23 f=[fvec1; y(1); y(2); y(3)]; 24
25 26 end
Sedan används båda dessa Matlabfunktioner för att lösa dierentialekva-tionerna med hjälp avode45. Detta görs genom kommandoraderna
1 clear 2 3 m_1=100; m_2=7.26; r=0.20; L1=1.70; L2=1.21; teta_4=pi/2; 4 X_0=L2*sin(teta_4); Z_0=L1+L2*cos(teta_4); 5 Ixx=((m_1)/12)*(3*r^2+4*L1^2); 6 Izz=(m_1*r^2)/2; 7 8 I_tot=[Ixx+m_2*Z_0^2 0 m_2*−X_0*Z_0; 0 ... Ixx+m_2*(X_0^2+Z_0^2) 0; m_2*−X_0*Z_0 0 Izz+m_2*X_0^2]; 9
10 y0=[0; 0; pi/24; 0; pi/6; 0]; T=3; 11 12 options=odeset('mass',@(t,y)lhs_hammer_anton(t,y,I_tot)); 13 14 [T,Y]=ode45(@(t,y)rhs_hammer_anton(t,y,I_tot),[0 ... T],y0,options); 15 16 plot(T,Y(:,4:6)) 17 18 legend('\theta_1(t)','\theta_2(t)','\theta_3(t)') 19 20 [m,n]=size(Y); 21 22 r0=[X_0 0 Z_0]'; 23 24 axis([−25 4 −4 4 0 4]) 25 hold on 26 dx = 0.15;xadd=0; 27 for k=1:m 28 for i=4:6 29 c(i−3)=cos(Y(k,i)); 30 s(i−3)=sin(Y(k,i)); 31 end 32 33 R=[c(1)*c(3)−c(2)*s(1)*s(3) −c(1)*s(3)−c(2)*c(3)*s(1) ... s(1)*s(2); c(3)*s(1)+c(1)*c(2)*s(3) ... c(1)*c(2)*c(3)−s(1)*s(3) c(1)*s(2); s(2)*s(3) ... c(3)*s(2) c(2)];
34 35 rk = R*r0; 36 37 xadd = xadd +dx; 38 rk(1) = rk(1)−xadd; 39 plot3(rk(1),rk(2),rk(3),'*') 40 41 pause(0.1) 42 end 43 hold off
Därplot(T,Y(:,4:6))ochlegendanvänds för att plotta ut vinklarna och
Bilaga C
Rörelseekvationer för
släggkastets rotationsfas
Här presenterar vi beräkningarna av Lagranges ekvationer, det vill säga rörel-seekvationerna för vår modell. Vi har valt att de niera tre stycken konstanter
k1 = l1l2m2, k2 = l22m2, k3 = l21m2
då de frekvent förekommer i beräkningarna. Ekvationen med avseende på q1 och ˙q1 ges av
d dt ∂L ∂ ˙q1 − ∂L ∂q1 = M1(t)
vilket efter beräkning ger
d dt ∂L ∂ ˙q1 −∂q∂L 1 = ¨q1 s 2 2(Ix0x0+ k3) + c22Iz0z0+ c23s22+ c22s24+ c24s22s23− 2c2c4s2s3s4 k2+ c4s22− c2s2s3s4 2k1) +¨q2 c3c24s2s3− c3s2s3− c2c3c4s4 k2− (c2c3s4) k1 +¨q3 c2Iz0z0+ c2s2 4− c4s2s3s4 k2− (s2s3s4) k1 + ¨q4((c3s2) k2+ (c3c4s2) k1) + ˙q22 c3c4s2s4+ c2c3c42s3− c2c3s3 k2+ (c3s2s4) k1 − ˙q2 3((c3s2s4) k1+ (c3c4s2s4) k2) − ˙q42(c3s2s4) k1 +2 ˙q1q˙2 c2s2(Ix0x0+ k3− Iz0z0) + s2 2s3s4+ 2c2c4s2− c22s3s4 k1+ c2c23s2− c2s2s24+ c2c24s2s23− c22c4s3s4+ c4s22s3s4 k2)
+2 ˙q1q˙3 c3c24s22s3− c3s22s3− c2c3c4s2s4 k2− c2c3s2s4k1 +2 ˙q1q˙4 c22c4s4− c4s22s23s4− c2c24s2s3+ c2s2s3s24 k2− s22s4+ c2c4s2s3 k1 + ˙q2q˙3 s2s23− c23s2+ c23c24s2− c24s2s23− s2s24 k2− s2Iz0z0 + ˙q2q˙4 c2c3s24− c2c3c24− c2c3− 2c3c4s2s3s4 k2 + ˙q3q˙4 s2s3s24+ 2c2c4s4− c42s2s3− s2s3 k2− 2c4s2s3k1 = M1(t).
Ekvationen med avseende på q2 och ˙q2 ges av
d dt ∂L ∂ ˙q2 − ∂L ∂q2 = M2(t)
vilket efter beräkning ger
d dt ∂L ∂ ˙q2 −∂q∂L 2 = ¨q1 c3c 2 4s2s3− c3s2s3− c2c3c4s4 k2− c2c3s4k1 +¨q2 Ix0x0 + k3+ c23c24+ s23 k2+ 2c4k1 −¨q3(c3s4k1+ c3c4s4k2) − ¨q4(c4s3k1+ s3k2) − ˙q2 1 c2s2(Ix0x0+ k3− Iz0z0) − c2 2s3s4− s22s3s4− 2c2c4s2 k1+ c2c24s2s23+ c2c23s2− c2s2s24− c22c4s3s4+ c4s22s3s4 k2) + ˙q32(c4s3s4k2+ s3s4k1) + ˙q42(s3s4k1) + ˙q1q˙3 s2Iz0z0+ 2c2c3s4k1+ s2s2 4+ s2s23− c23s2− c23c24s2+ 2c2c4s3s4− c24s2s23 k2 + ˙q1q˙4 c2c3s24− c2c3c24− c2c3− 2c3c4s3s4 k2− 2c2c3c4k1 +2 ˙q2q˙3 c3s3− c3c24s3 k2 − 2 ˙q2q˙4 s4k1+ c23c4s4k2 − ˙q3q˙4 2c3c4k1+ c3c24− c3s24+ c3 k2 + s2 m1gl21 + m2gl1 = M2(t).
Ekvationen med avseende på q3 och ˙q3 ges av
d dt ∂L ∂ ˙q3 − ∂L ∂q3 = M3(t)
vilket efter beräkning ger d dt ∂L ∂ ˙q3 − ∂L ∂q3 = ¨q1 c2Iz 0z0+ c2s2 4− c4s2s3s4 k2− s2s3s4k1 −¨q2(c3c4k1+ c3c4s4k2) + ¨q3 Iz0z0 + s2 4k2 − ˙q2 1 c3c24s22s3− c3s22s3− c2c3c4s2s4 k2− (c2c3s2s4) k1 − ˙q22 c3s3− c3c24s3 k2 − ˙q1q˙2 s2Iz0z0+ s2s24+ c24s2s23+ c23s2− c23c24s2− s2s23− 2c2c4s3s4 k2− 2c2s3s4k1 + ˙q1q˙4 s2s3s24+ s2s3− c24s2s3+ 2c2c4s4 k2+ ˙q2q˙4 c4s24+ c3− c3c24 k2 + ˙q3q˙4(2c4s4k2) = M3(t).
Ekvationen med avseende på q4 och ˙q4 ges av
d dt ∂L ∂ ˙q4 − ∂L ∂q4 = M4(t)
vilket efter beräkning ger
d dt ∂L ∂ ˙q4 − ∂q∂L 4 = ¨q1(c3c4s2k1+ c3s2k2) −¨q2(c4s3k1+ s3k2) + ¨q4k2 − ˙q2 1 c22c4s4− c4s22s23s4+ c2s2s3s42− c2c24s2s3 k2− c2c4s2s3+ s22s4 k1 + ˙q22 s4k1+ c23c4s4k2 − ˙q23(c4s4k2) + ˙q1q˙2 c2c3+ c2c3c24− c2c3s24+ 2c3c4s2s3s4 k2+ 2c2c3c4k1 − ˙q1q˙3 s2s3+ s2s3s42− s2s3c24+ 2c2c4s4 k2 − ˙q2q˙3 c3+ c3s24− c3c24 k2+ s4(m2gl2) = M4(t).
