Hemelmechanica


(Credits voor de foto’s van de hemellichamen: NASA)

Newton

In het heelal is het zo dat alles beweegt ten opzichte van elkaar dus we gaan eerst wat referenties definiëren. Bovendien is het zo dat in het heelal de zwaartekracht regeert. Mijn uitgangspunt is dan ook de zwaartekrachtwet van Newton:

Hierin zijn m1 en m2 twee massa’s, r is de afstand tussen (de zwaartepunten van) die massa’s en G is de gravitatieconstante:
Om te beginnen vormen we een coördinatenstelsel x-y-z waarvan de assen loodrecht op elkaar staan. De oorsprong van dit stelsel valt samen met het zwaartepunt, zeg maar het middelpunt, van de Zon. Zon-gerelateerde parameters geef ik de index 0 mee, net alsof de Zon het nulde object, of de nulde planeet, is. De coördinaten x0, y0 en z0 geven aan waar de Zon zich bevindt en die zijn dus alledrie nul omdat we net hebben afgesproken dat het zwaartepunt van de Zon samenvalt met de oorsprong. Op een afstand r10 (de afstand van object 1 tot object 0) van de oorsprong bevindt zich een planeet p1 met massa m1 die in een baan om de Zon draait. Deze planeet bevindt zich op coördinaten (x1, y1, z1). Er geldt:
De zwaartekracht die de planeet op de Zon uitoefent is (opgesplitst in de x, y en z richting):


Een tweede planeet p2 met massa m2 op afstand r20 van de oorsprong oefent op dezelfde wijze zwaartekracht uit op de Zon:


Indien er n planeten zijn dan is de totale zwaartekrachtwerking van deze planeten op de Zon:


Door al dit getrek aan de Zon beweegt de Zon uiteraard, maar niet binnen het x-y-z-stelsel want de oorsprong van dit stelsel hebben we immers in het zwaartepunt van de Zon neergelegd. Vervolgens stellen we ons voor dat er zich ergens een ‘vast punt’ in de ruimte bevindt en dat punt bevindt zich in de oorsprong van het coördinatenstelsel ξ-η-ζ. Ten opzichte van het ξ-η-ζ-stelsel beweegt de Zon natuurlijk wel, en wel onder invloed van de hierboven beschreven krachten:


Door de vergelijkingen (5) en (6) te combineren ontstaat:


Omdat zwaartekracht betekent dat alles aan alles trekt, trekt p2 niet alleen aan de Zon maar ook aan p1:


En dat geldt ook voor de planeten p3, p4, ..., pn:


Nog even voor de duidelijkheid:
En:
Bovendien trekt de Zon ook nog aan p1:


Alles wat er aan p1 trekt is dus:


Hierbij dienen we steeds te bedenken dat x0, y0 en z0 alledrie nul zijn omdat de Zon zich in de oorsprong bevindt. Daarnaast valt de term waarbij i = 1 er vanzelf uit omdat de teller van de breuk dan nul wordt. De bovenstaande vergelijkingen geven de krachten die p1 laten bewegen en die kan ik ook schrijven als:


Vanuit het ξ-η-ζ-stelsel bezien wordt dit:


Voor de tweede afgeleiden van de ξ-η-ζ-coördinaten kan ik, met behulp van de vergelijkingen (7), ook schrijven:


Nu grijpen we weer eens in de trucendoos en we stellen:
Oftewel:
De functie Θ ga ik differentiëren naar respectievelijk x1, y1 en z1:


Met de bovenstaande afgeleiden kan ik de vergelijkingen (16) ook als volgt opschrijven:


Ik stel:

Zo komen we tot deze set vergelijkingen:


De functie Θ heeft alle verstorende invloeden, door de andere planeten op p1, in zich en daarom noemen we Θ de verstoringsfunctie. Het is ook belangrijk om op te merken dat ik de planeten weliswaar genummerd heb als p1, p2, ..., pn, maar dat dit volkomen los staat van hoe de werkelijke volgorde van de planeten is. De vergelijkingen (22) zijn de bewegingsvergelijkingen van planeet p1, inclusief alle verstorende invloeden, maar ik kan dus iedere planeet aanmerken als p1. p1 is niet noodzakelijkerwijs de binnenste planeet, maar een willekeurige planeet! Daarom noem ik x1 vanaf nu x, y1 wordt y, z1 wordt z en r10 wordt r. De vergelijkingen (22) worden dan:


Merk op dat wanneer er geen andere planeten dan p1 aanwezig zijn, dat p1 dan ongestoord zijn rondjes kan draaien en de vergelijkingen (23) overgaan in:


Door de vergelijking (24a) met y te vermenigvuldigen en (24b) met x, en vervolgens de resultaten van elkaar af te trekken krijgen we dit:
Ditzelfde ga ik ook doen door (24a) met z te vermenigvuldigen en (24c) met x, en (24b) met z te vermenigvuldigen en (24c) met y:

De vergelijkingen (25) ga ik integreren:


Hierin zijn de c’s integratieconstanten. De vergelijkingen (26) ga ik vermenigvuldigen met respectievelijk z, −y en x, en ik tel ze alledrie bij elkaar op:
Dit is de vergelijking van een plat vlak, dus de planeet beweegt in een plat vlak en de Zon bevindt zich in datzelfde vlak (in de oorsprong). Dit resultaat gaan we straks gebruiken. Nog even dit tussendoortje, voor r2 geldt:
Hieruit volgt door te differentiëren:
Hier ga ik gebruik van maken bij het volgende. Ik ga de vergelijkingen (24) vermenigvuldigen met respectievelijk 2dx, 2dy en 2dz, en bij elkaar optellen:
Vergelijking (30) ga ik vervolgens integreren:
Dit resultaat parkeren we even. Nu ga ik de vergelijkingen (26) kwadrateren en bij elkaar optellen:
Ik stel:
Waarmee vergelijking (32) tenslotte wordt:

Pythagoras

Indien de planeet een infinitesimaal stukje aflegt van zijn baan, dan geldt er voor dat stukje baanlengte in Cartesische coördinaten (gewone x-y-z-coördinaten):

Dat is simpelweg de stelling van Pythagoras toepassen.


Figuur 1

Maar in poolcoördinaten ligt het iets ingewikkelder:


Figuur 2

Laat ik even inzoomen:

De infinitesimale verandering van r is dr, en loodrecht daarop staat r dθ. Het stukje baanlengte in poolcoördinaten wordt dan:
Uit (35) en (36) volgt:
Dit vul ik in in (34):
Of iets anders geschreven:
Nu ga ik vergelijking (31) iets anders opschrijven:
En vergelijking (32) ook:
Door (40) en (41) te combineren krijg ik een uitdrukking voor dt:
En die vul ik in in (38):
Door dit iets anders op te schrijven heb ik in één klap de afgeleide van r (naar θ):
En door (44) nul te stellen vind ik de extreme waarden van r:
De oplossing r = 0 is die waarbij beide hemellichamen op ramkoers liggen en waarbij ze uiteindelijk op elkaar botsen. Maar het zijn natuurlijk de andere oplossingen die ons interesseren. Ik stel:

Waarmee ik (45) als volgt kan schrijven:
De punten Q en q zijn de apsiden, de punten van maximale en minimale afstand tot de Zon. Q is de afstand van de Zon tot het aphelium (het punt van grootste afstand tot de Zon):
En q is de afstand van de Zon tot het perihelium (het punt van kleinste afstand tot de Zon):
Met deze kennis ga ik vergelijking (43) omschrijven:
Vervolgens ga ik integreren:
Hier heb ik de vergelijking van een ellips! De integratieconstante θ0 is de hoek van waaraf θ wordt bepaald en die stel ik voor het gemak gelijk aan nul. Daardoor kan ik vergelijking (52) tenslotte schrijven als (waarbij ik gebruik maak van de vergelijkingen (46) en (47) om uiteindelijk c4, C en k kwijt te raken):
De hoek θ is de werkelijke anomalie. We kennen ε als de lineaire excentriciteit van een ellips en deze loopt van nul (in het geval van een perfecte cirkel: a = b) tot a (wanneer de ellips zo plat is als een dubbeltje: b = 0):
Daarnaast kennen we e als de numerieke excentriciteit, of kortweg excentriciteit, van een ellips, deze loopt van nul (in het geval van een perfecte cirkel: a = b) tot één (wanneer de ellips zo plat is als een dubbeltje: b = 0):
De halve lange as van de ellips is a, en b is de halve korte as. De teller van (53) kan ik dus omschrijven als volgt:
Ik stel:
Waarmee (53) tenslotte wordt:
Ik schrijf nu (47) even iets anders op en vind zo voor C:
En vervolgens ga ik verder met vergelijking (42) en ik vul hier de waarden van de integratieconstanten in volgens (46) en (59), en ik maak ook gebruik van (21):
Bovenstaande vergelijking ga ik integreren waarbij ik als integratiegrenzen r = a (1 − e) en r = a (1 + e) neem, oftewel de punten van kleinste afstand (het perihelium) en grootste afstand (het aphelium) tot de Zon. De tijdsduur hiervan (om van perihelium tot aphelium te komen) is uiteraard de helft van de totale omlooptijd T:
En zo hebben we de omlooptijd T gevonden. Nu willen we nog de snelheid weten. Voor de snelheid v geldt in het algemeen:
Vanuit vergelijking (40) kan ik dus een uitdrukking vinden voor de snelheid:
Het is wel interessant (en handig, dat zullen we later merken) om deze snelheid te ontbinden in een radiële (naar buiten gericht, in de richting van r) component vr en een laterale (met de draairichting, loodrecht op r) component vl. De radiële component vind ik door r te differentiëren naar de tijd:
De afgeleide van θ naar de tijd haal ik uit vergelijking (39) en vervolgens maak ik ook nog gebruik van (58) en (59):
De laterale snelheid en de radiële snelheid vormen samen uiteraard de totale snelheid (middels de stelling van Pythagoras):
Hiermee vind ik voor de laterale snelheid:
De hoek tussen de raaklijn aan de ellips (de richting van de snelheid) en het verlengde van r (de richting van de radiële snelheid) noem ik ψ:

Figuur 3
Ik vind ψ als volgt:

Hieruit volgt voor de sinus respectievelijk de cosinus van ψ:

Ik noem ook nog de hoeksnelheid, dit is de laterale snelheid gedeeld door r:
Vergelijking (38) heeft ook nog iets interessants in petto:
Het linkerlid is een infinitesimaal oppervlakje dO (de oppervlakte van een gelijkbenige driehoek is halve basis maal hoogte, zie ook figuur 1) van de ellips en het rechterlid is de tijd voorafgegaan door een aantal constanten. Door dit te gaan integreren tussen willekeurige grenzen volgt hieruit:

Kepler

Hier staat dat een planeet in gelijke tijdsintervallen gelijke oppervlakten bestrijkt. In het Engels noemen ze dit de “areal velocity”, hetgeen letterlijk vertaald oppervlaktesnelheid betekent maar eigenlijk hebben we hier in het Nederlands geen woord voor. Wij kennen dit alleen als de perkenwet van Kepler. Ik ga hier gewoon oppervlaktesnelheid introduceren:

Hiernaast de originele afleiding van Newton van de perkenwet van Kepler. Zoals bekend was Newton (ook) heel goed in meetkunde, waarvan akte.

Aldus hebben we de ‘gewone’ (of lineaire) snelheid v (vergelijking (63)), en de oppervlaktesnelheid v-met-een-streepje-erboven (vergelijking (73)). Die hebben we nog nodig om de integratieconstanten c1, c2 en c3 te vinden. Na al het voorgaande wiskundige werk heb ik nu eerst even een plaatje gemaakt voor onze planeet, de Aarde:

Figuur 4
Zoals gezegd gebruik ik telkens een nul als index voor Zon-gerelateerde variabelen, een één voor ‘de-eerste-planeet’-gerelateerde variabelen en een twee (of hoger) voor de overige planeten. Indien a = b dan hebben we een cirkel en zijn e en ε beide nul. En indien b = 0 dan is de ellips zo plat als een dubbeltje en is e = 1 en ε = a.

De twee brandpunten van een ellips liggen op een afstand ε (= ae) van het middelpunt van een ellips. Indien ik een lijn trek vanaf een brandpunt naar de top van de ellips dan is de lengte van die lijn:
Zoals je kunt zien in het onderstaande plaatje:

Figuur 5
Dit is ook helemaal consistent met de definitie van een ellips, want dat is “de verzameling der punten, waarvoor de som der afstanden tot de beide brandpunten constant is (en wel gelijk aan de lange as van de ellips = 2a)”. Ik heb ook een hoek φ aangegeven, de excentriciteitshoek, de sinus hiervan is:
De sinus van de complementaire hoek van φ is dan (of de cosinus van φ, dat komt op het zelfde neer):
Hierdoor krijgen we de beschikking over enkele handige relaties tussen de excentriciteit en de complementaire excentriciteit:



Maar hoe geef je dit alles een positie in de ruimte? Om te beginnen nemen we het vlak waar de baan van de Aarde in ligt, dit noemen we het eclipticavlak. De Aarde en de Zon bevinden zich dus per definitie in het eclipticavlak. Vanaf de Aarde gezien zie je de Zon door het eclipticavlak trekken en vanaf de Zon gezien zie je de Aarde door het eclipticavlak trekken. Het eclipticavlak is het referentievlak van ons zonnestelsel. Waar de ‘randen’ van het eclipticavlak de hemelbol ‘raken’ noemen we de ecliptica, dat is dus een denkbeeldige cirkel op de hemelbol.

De evenaar van de Aarde noemen we met een mooi woord ook wel equator. De projectie van de equator op de hemelbol is de hemelequator. Wanneer je op een willekeurige plaats op de evenaar staat en je kijkt recht omhoog (naar het zenit, dat is per definitie het punt recht boven je, en nadir is per definitie het punt recht onder je) dan kijk je dus altijd tegen de hemelequator aan (maar die zie je natuurlijk niet, want ook de hemelequator is een puur denkbeeldige cirkel). De aardas staat niet loodrecht op het eclipticavlak maar helt ongeveer 23.5 graden over. Het vlak waar de evenaar en de hemelequator in liggen, het equatorvlak, staat dus onder een hoek van 23.5 graden met het eclipticavlak. De ecliptica en de hemelequator maken uiteraard dezelfde hoek met elkaar want ze liggen in het eclipticavlak respectievelijk het equatorvlak.

Er is ook nog een derde vlak, het baanvlak. Dit is het vlak waar de baan van een planeet in ligt. Voor de Aarde vallen het baanvlak en het eclipticavlak per definitie samen. De baanvlakken van de andere planeten (in ons zonnestelsel) maken allemaal een (kleine) hoek met het eclipticavlak, deze hoek heet de inclinatiehoek i (voor de Aarde geldt altijd i = 0), of kortweg inclinatie of baanhelling. Van opzij gezien ziet dat er zo uit:

Figuur 6
Het eclipticavlak stellen we als x-y-vlak en vergelijking (27) liet ons zien dat het baanvlak van een planeet beschreven wordt als:
Waaruit na enig denkwerk volgt dat c1 tweemaal de projectie is van de oppervlaktesnelheid:
Een planeetbaan gaat op twee plaatsen door het eclipticavlak heen en deze twee plaatsen noemen we knopen (het snijpunt van de blauwe lijn en de groene lijn in figuur 6). Wanneer de planeet aan het ‘klimmen’ is terwijl die door het eclipticavlak heen gaat dan spreken we van de klimmende knoop en wanneer de planeet aan het ‘dalen’ is dan spreken we van, je raadt het al, de dalende knoop. De hoek tussen de klimmende knoop en het perihelium heet ω. Deze hoek noemen we het argument van het perihelium. Van bovenaf gezien (we kijken dan op het baanvlak) ziet dat er zo uit (waarbij ik de knopen op willekeurige plaatsen heb getekend, want dat verschilt natuurlijk per planeetbaan):

Figuur 7
De blauwe lijn is de lijn door de beide knopen, de knopenlijn, en dit is tevens de snijlijn van baanvlak en eclipticavlak. Het deel van het baanvlak rechts van de knopenlijn ligt boven het eclipticavlak en het deel van het baanvlak links van de knopenlijn ligt onder het eclipticavlak (ervan uitgaande dat de planeet tegen de wijzers van de klok in draait). Als een planeet ‘de andere kant op draait’, dat heet retrograde, dan wisselen in het plaatje hierboven de klimmende knoop en de dalende knoop van plaats. Een 3D-impressie ziet er als volgt uit.

Figuur 8
Na nog meer denkwerk volgt dat c2 en c3 ook tweemaal de projectie zijn van de oppervlaktesnelheid:

Waardoor vergelijking (27), de vergelijking van het baanvlak, tenslotte wordt:
Let wel, dit is de vergelijking van het baanvlak en niet van de baan zelf. Indien de inclinatie van het baanvlak nul graden is (i = 0), dan is de cosinus van i gelijk aan één en de sinus van i gelijk aan nul en reduceert (81) inderdaad tot z = 0. En stel dat de inclinatie van het baanvlak negentig graden is (i = π/2) dan staat het baanvlak loodrecht op het eclipticavlak en wordt de vergelijking van het baanvlak y = 0 voor ω is nul graden en x = 0 voor ω is negentig graden.

Op de hemelbol hebben we dus twee denkbeeldige cirkels, de ecliptica en de hemelequator. Vanaf de Aarde gezien zien we de Zon de ecliptica volgen. Deze twee cirkels hebben noodzakerlijkerwijs twee snijpunten, namelijk daar waar de Zon recht boven de evenaar staat. Die twee snijpunten heten het lentepunt, de Zon is dan aan het ‘klimmen’, en het herfstpunt, de Zon is dan aan het ‘dalen’. Men heeft voor het lentepunt gekozen om uiteindelijk de planeetbaan ‘vast’ te leggen in de ruimte middels de hoek Ω. Deze hoek noemen we de lengte van de klimmende knoop.

Figuur 9
Het blauwe vlak is het eclipticavlak (maar dat had ik net zo goed rechthoekig kunnen tekenen). De zwarte lijn, de knopenlijn, is de snijlijn van baanvlak en eclipticavlak. De groene lijn is de snijlijn van equatorvlak en eclipticavlak. Met de volgende parameters, de baanelementen, kan ik aldus een planeetbaan (binnen ons zonnestelsel) volledig vastleggen: En de zesde parameter, de epoche, hoort er wel degelijk bij, want omdat alles aan alles ‘loopt te trekken’ (zwaartekracht) hebben de baanelementen geenszins vaste waarden maar variëren ze (langzaam) in de tijd. Bij de integraal van vergelijking (61) heb ik als ondergrens t = 0 genomen, maar wat betekent t = 0? De epoche is dat moment waarvoor t = 0 en waarvoor de overige baanelementen gelden en wordt uitgedrukt als datum, bijvoorbeeld 1 juli 2015. Voor de actuele positie van een planeet hebben we ook nog de hoek θ nodig en de afstand r (zie figuur 4). Indien θ = 0 staat de planeet in het perihelium (r = minimaal) en als θ = π staat de planeet in het aphelium (r = maximaal).

Ik ga dit onderdeel hemelmechanica compleet maken met de volgende parameters die ook veel gebruikt worden. Allereerst de gemiddelde dagelijkse beweging van een planeet:
Doorgaans wordt n uitgedrukt in graden per dag. Je kunt n dus zien als een soort gemiddelde hoeksnelheid, maar niet echt natuurlijk want de werkelijke gemiddelde hoeksnelheid is:
De hoek θ is de werkelijke anomalie, maar daarnaast kennen we ook de excentrische anomalie E. Ik heb de hoek E aangegeven in onderstaande figuur:

Figuur 10
Om de ellips heb ik een cirkel getekend met straal a. Door het verschil a − r (het gele gebied is de verzameling van alle (a − r)’s) af te meten, langs r, vanaf het middelpunt M kom ik bij het punt K. De afstand MK is dus a − r. De cosinus van E vind ik als volgt:
Waaruit volgt:
Of, met behulp van vergelijking (58):
Of, door (85a) te kwadrateren en daar het kwadraat van (85b) van af te trekken (met hulp van 77b):
Tenslotte ga ik nogmaals de integraal van vergelijking (61) uitvoeren, maar nu als onbepaalde integraal:
Vervolgens maak ik gebruik van (82) en (85). En (61) heeft ons reeds laten zien dat c5 = 0:
Het linkerlid noemen we de middelbare anomalie M:
Waardoor we (87) dus kunnen schrijven als:
Verder wil ik nog opmerken dat in tabellen met baanelementen regelmatig in plaats van het argument van het perihelium de lengte van het perihelium gegeven wordt: ϖ. Deze zijn in elkaar om te rekenen als volgt:

Met de kennis die we nu hebben kan ik (39) nog iets anders opschrijven, dat is ook wel handig voor later gebruik:
Het impulsmoment is ook goed om paraat te hebben, dat is (in dit geval) het product van r met de snelheid die daar loodrecht op staat (de laterale snelheid dus). Die kennen we binnen de astronomie in twee smaken, allereerst de massaloze variant:
En het echte impulsmoment, dus met massa:

Dan hebben we nu een stevige basis hemelmechanica waarbij we ervan uitgegaan zijn dat er geen storende invloeden zijn. Voor een planeet bestaan er dus zes baanelementen: a (halve lange baanas), e (excentriciteit), i (inclinatie), ω (argument van het perihelium), Ω (lengte van de klimmende knoop) en t (epoche). Hiermee is zo ongeveer alles uit te rekenen wat je van een planeet zou willen weten. Allereerst de constante k:
Voor p geldt:
De afstand r als functie van de werkelijke anomalie θ:
De afstand r als functie van de excentrische anomalie E:
In combinatie met de werkelijke anomalie θ:

De afstand van de Zon tot het aphelium:
De afstand van de Zon tot het perihelium:
De (numerieke) excentriciteit:
De complementaire excentriciteit:
De lineaire excentriciteit:
De omlooptijd:
De snelheid:
De radiële snelheid:
De laterale snelheid:
De hoeksnelheid:
De oppervlaktesnelheid:
De hoek ψ, tussen het verlengde van r en de richting van de snelheid:

De excentriciteitshoek φ:

De vergelijking van het baanvlak:
De gemiddelde dagelijkse beweging:
De excentrische anomalie:
De middelbare anomalie:
De lengte van het perihelium:
Het impulsmoment, zonder massa en met massa:

En dit alles volgt uit één simpel uitgangspunt: de zwaartekrachtwet van Newton!