Wie man Neptun aus der Uranus-Umlaufbahn entdeckt (durch Computersimulation)

11

Ich möchte die Existenz eines anderen Planeten (Neptun) demonstrieren, indem ich die Diskrepanz zwischen der Beobachtung der Uranus-Umlaufbahn und der mathematischen Vorhersage untersuche. Diese Arbeit wurde von Le Verrier gemacht und ich möchte seine Methode verstehen.

Ich habe das Kapitel 2 "Die Entdeckung des Neptun (1845-1846)" in der Biographie "Le Verrier - Prächtiger und verabscheuungswürdiger Astronom" gelesen, aber es ist zu ausführlich und ich habe seine Arbeit nicht sehr gut verstanden.

Ich untersuche das Drei-Körper-Problem (Sonne, Uranus, Neptun) über Matlab und das Zwei-Körper-Problem (Sonne, Uranus), wobei ich den Ausgangszustand von hier aus nehme:

http://nssdc.gsfc.nasa.gov/planetary/factsheet/uranusfact.html

Ich habe diese Methode ausprobiert: Ich habe Uranus mit dem Max in das Perihel gelegt. Die Orbitalgeschwindigkeit und ich berechnen die Semi-Major-Achse und sie ist genauer als diejenige, die erhalten wird, wenn Uranus und Neptun mit ihren jeweiligen max. Umlaufgeschwindigkeit.

Hier ein cooles Bild mit Matlab: Hier ein cooles Bild

Kann mir jemand helfen? Was muss ich tun und mit welchen Daten muss ich meine Vorhersage vergleichen? Auch ein einfacher Link könnte hilfreich sein.

Sergio Piccione
quelle

Antworten:

11

Folgendes habe ich getan:

  • Aufgrund ihrer Masse ist es am sichersten, zunächst Jupiter und Saturn sowie Uranus zu betrachten. Es könnte auch fruchtbar sein, die Erde in die Analyse einzubeziehen, relative Positionen, Beobachtungswinkel usw. zu erhalten. Also werde ich überlegen:
    • Sonne
    • Erde
    • Jupiter
    • Saturn
    • Uranus
    • Neptun
  • Holen Sie sich die Standard-Gravitationsparameter (μ) für alle
  • Erhalten Sie Anfangspositionen und Geschwindigkeiten über JPL / HORIZONS für alle diese Planeten. Ich hatte einige Daten von J2000.5 herumliegen, also habe ich nur die Zustandsvektoren vom 1. Januar 2000 um 12 Uhr verwendet.
  • Schreiben Sie einen N-Body-Integrator mit integrierten MATLAB-Tools. Integrieren Sie dieses unvollständige Sonnensystem einmal ohne Neptun und einmal mit Neptun.
  • Analysieren und vergleichen!

Also, hier sind meine Daten und der N-Körper-Integrator:

function [t, yout_noNeptune, yout_withNeptune] = discover_Neptune()

    % Time of integration (in years)
    tspan = [0 97] * 365.25 * 86400;

    % std. gravitational parameters [km/s²/kg]
    mus_noNeptune = [1.32712439940e11; % Sun
                     398600.4415       % Earth
                     1.26686534e8      % Jupiter
                     3.7931187e7       % Saturn
                     5.793939e6];      % Uranus

    mus_withNeptune = [mus_noNeptune
                       6.836529e6]; % Neptune

    % Initial positions [km] and velocities [km/s] on 2000/Jan/1, 00:00
    % These positions describe the barycenter of the associated system,
    % e.g., sJupiter equals the statevector of the Jovian system barycenter.
    % Coordinates are expressed in ICRF, Solar system barycenter
    sSun     = [0 0 0 0 0 0].';
    sEarth   = [-2.519628815461580E+07  1.449304809540383E+08 -6.175201582312584E+02,...
                -2.984033716426881E+01 -5.204660244783900E+00  6.043671763866776E-05].';
    sJupiter = [ 5.989286428194381E+08  4.390950273441353E+08 -1.523283183395675E+07,...
                -7.900977458946710E+00  1.116263478937066E+01  1.306377465321731E-01].';
    sSaturn  = [ 9.587405702749230E+08  9.825345942920649E+08 -5.522129405702555E+07,...
                -7.429660072417541E+00  6.738335806405299E+00  1.781138895399632E-01].';
    sUranus  = [ 2.158728913593440E+09 -2.054869688179662E+09 -3.562250313222718E+07,...
                 4.637622471852293E+00  4.627114800383241E+00 -4.290473194118749E-02].';
    sNeptune = [ 2.514787652167830E+09 -3.738894534538290E+09  1.904284739289832E+07,...
                 4.466005624145428E+00  3.075618250100339E+00 -1.666451179600835E-01].';

    y0_noNeptune   = [sSun; sEarth; sJupiter; sSaturn; sUranus];
    y0_withNeptune = [y0_noNeptune; sNeptune];

    % Integrate the partial Solar system 
    % once with Neptune, and once without
    options = odeset('AbsTol', 1e-8,...
                     'RelTol', 1e-10);

    [t, yout_noNeptune]   = ode113(@(t,y) odefcn(t,y,mus_noNeptune)  , tspan, y0_noNeptune  , options);
    [~, yout_withNeptune] = ode113(@(t,y) odefcn(t,y,mus_withNeptune),     t, y0_withNeptune, options);

end

% The differential equation 
%
%    dy/dt = d/dt [r₀ v₀ r₁ v₁ r₂ v₂ ... rₙ vₙ]    
%          = [v₀ a₀ v₁ a₁ v₂ a₂ ... vₙ aₙ]    
%
%  with 
%
%    aₓ = Σₘ -G·mₘ/|rₘ-rₓ|² · (rₘ-rₓ) / |rₘ-rₓ| 
%       = Σₘ -μₘ·(rₘ-rₓ)/|rₘ-rₓ|³  
%
function dydt = odefcn(~, y, mus)

    % Split up position and velocity
    rs = y([1:6:end; 2:6:end; 3:6:end]);
    vs = y([4:6:end; 5:6:end; 6:6:end]);

     % Number of celestial bodies
    N = size(rs,2);

    % Compute interplanetary distances to the power -3/2
    df  = bsxfun(@minus, permute(rs, [1 3 2]), rs);
    D32 = permute(sum(df.^2), [3 2 1]).^(-3/2);
    D32(1:N+1:end) = 0; % (remove infs)

    % Compute all accelerations     
    as = -bsxfun(@times, mus.', D32);              % (magnitudes)    
    as = bsxfun(@times, df, permute(as, [3 2 1])); % (directions)    
    as = reshape(sum(as,2), [],1);                 % (total)

    % Output derivatives of the state vectors
    dydt = y;
    dydt([1:6:end; 2:6:end; 3:6:end]) = vs;
    dydt([4:6:end; 5:6:end; 6:6:end]) = as;

end

Hier ist das Treiberskript, mit dem ich einige nette Handlungen veröffentlicht habe:

clc
close all

% Get coordinates from N-body simulation
[t, yout_noNeptune, yout_withNeptune] = discover_Neptune();

% For plot titles etc.
bodies = {'Sun'
          'Earth'
          'Jupiter'
          'Saturn'
          'Uranus'
          'Neptune'};


% Extract positions
rs_noNeptune   = yout_noNeptune  (:, [1:6:end; 2:6:end; 3:6:end]);
rs_withNeptune = yout_withNeptune(:, [1:6:end; 2:6:end; 3:6:end]);



% Figure of the whole Solar sysetm, just to check
% whether everything went OK
figure, clf, hold on
for ii = 1:numel(bodies)
    plot3(rs_withNeptune(:,3*(ii-1)+1),...
          rs_withNeptune(:,3*(ii-1)+2),...
          rs_withNeptune(:,3*(ii-1)+3),...
          'color', rand(1,3));
end

axis equal
legend(bodies);
xlabel('X [km]');
ylabel('Y [km]');
title('Just the Solar system, nothing to see here');


% Compare positions of Uranus with and without Neptune
rs_Uranus_noNeptune   = rs_noNeptune  (:, 13:15);
rs_Uranus_withNeptune = rs_withNeptune(:, 13:15);

figure, clf, hold on

plot3(rs_Uranus_noNeptune(:,1),...
      rs_Uranus_noNeptune(:,2),...
      rs_Uranus_noNeptune(:,3),...
      'b.');

plot3(rs_Uranus_withNeptune(:,1),...
      rs_Uranus_withNeptune(:,2),...
      rs_Uranus_withNeptune(:,3),...
      'r.');

axis equal
xlabel('X [km]');
ylabel('Y [km]');
legend('Uranus, no Neptune',...
       'Uranus, with Neptune');


% Norm of the difference over time
figure, clf, hold on

rescaled_t = t/365.25/86400;

dx = sqrt(sum((rs_Uranus_noNeptune - rs_Uranus_withNeptune).^2,2));
plot(rescaled_t,dx);
xlabel('Time [years]');
ylabel('Absolute offset [km]');
title({'Euclidian distance between'
       'the two Uranuses'});


% Angles from Earth
figure, clf, hold on

rs_Earth_noNeptune   = rs_noNeptune  (:, 4:6);
rs_Earth_withNeptune = rs_withNeptune(:, 4:6);

v0 = rs_Uranus_noNeptune   - rs_Earth_noNeptune;
v1 = rs_Uranus_withNeptune - rs_Earth_withNeptune;

nv0 = sqrt(sum(v0.^2,2));
nv1 = sqrt(sum(v1.^2,2));

dPhi = 180/pi * 3600 * acos(min(1,max(0, sum(v0.*v1,2) ./ (nv0.*nv1) )));
plot(rescaled_t, dPhi);

xlabel('Time [years]');
ylabel('Separation [arcsec]')
title({'Angular separation between the two'
       'Uranuses when observed from Earth'});

was ich hier Schritt für Schritt beschreiben werde.

Zunächst ein Diagramm des Sonnensystems, um zu überprüfen, ob der N-Körper-Integrator ordnungsgemäß funktioniert:

das Sonnensystem

Nett! Als nächstes wollte ich den Unterschied zwischen den Positionen von Uranus mit und ohne den Einfluss von Neptun sehen. Also habe ich nur die Positionen dieser beiden Uranus extrahiert und sie aufgezeichnet:

Zwei Uranus mit und ohne Neptun

... das ist kaum sinnvoll. Selbst wenn Sie stark hineinzoomen und das Heck herausdrehen, ist dies keine nützliche Handlung. Also habe ich mir die Entwicklung der absoluten euklidischen Distanz zwischen den beiden Uranussen angesehen:

Zeitliche Entwicklung der euklidischen Distanz zwischen den beiden Uranus

Das sieht jetzt eher so aus! Etwa 80 Jahre nach Beginn unserer Analyse sind die beiden Uranus fast 6 Millionen km voneinander entfernt!

So groß das auch klingen mag, im größeren Maßstab der Dinge könnte dies im Lärm ertrinken, wenn wir hier auf der Erde Messungen durchführen. Außerdem erzählt es immer noch nicht die ganze Geschichte, wie wir gleich sehen werden. Schauen wir uns als nächstes den Winkeldifferenz zwischen den Beobachtungsvektoren von der Erde zu den beiden Uranus an, um zu sehen, wie groß dieser Winkel ist und ob er über den Beobachtungsfehlerschwellen hervorstechen kann:

Winkeltrennung zwischen den beiden Uranus

... whoa! Weit über 300 Bogensekunden Unterschied, plus alle Arten von wackeligen Bobbley Timey Wimey Ripple. Das scheint gut in den Beobachtungsmöglichkeiten der Zeit zu liegen (obwohl ich so schnell keine verlässliche Quelle dafür finden kann; irgendjemand?)

Nur zum guten Teil habe ich auch diese letzte Handlung produziert, die Jupiter und Saturn aus dem Bild lässt. Obwohl einige Störungstheorie in der 17 entwickelt worden war th und 18 th Jahrhundert wurde es nicht sehr gut entwickelt und ich bezweifle sogar Le Verrier Jupiter in Betracht zog (aber auch hier könnte ich falsch sein, bitte korrigieren Sie mich , wenn Sie mehr wissen).

Also, hier ist die letzte Handlung ohne Jupiter und Saturn:

Winkeltrennung zwischen den beiden Uranus, so dass Jupiter und Saturn aus der Gleichung herausfallen

Obwohl es Unterschiede gibt, sind sie winzig und vor allem für die Entdeckung von Neptun irrelevant.

Rody Oldenhuis
quelle
Geniale Antwort!
Zephyr
4

Wenn ich das richtig verstehe, modellieren Sie die Umlaufbahn von Uranus als Ellipse und möchten sie mit der tatsächlichen Umlaufbahn von Uranus vergleichen, die von Neptun gestört wird? Ich habe keine Antwort, aber wo kann ich Planeten / Sterne / Monde / usw. Positionen finden / visualisieren? erklärt, wie SPICE, HORIZONS und andere Tools verwendet werden, um die wahre Position von Uranus zu einem bestimmten Zeitpunkt in + -15000 Jahren zu ermitteln, einschließlich der am besten passenden elliptischen Parameter (mithilfe der HORIZONS-Funktion "Orbitalelemente").

Natürlich wird alles, was Sie tun, in gewissem Sinne "kreisförmig" sein, da HORIZONS berechnete Position von Uranus in der Vergangenheit bereits Neptuns Störungen enthält.

Wenn Sie Tabellen mit Uranus-Positionsvorhersagen oder etwas aus der Vergangenheit finden könnten, hätten Sie möglicherweise etwas.

Übrigens, zögern Sie nicht, mich zu kontaktieren (siehe Profil für Details), wenn dieses Projekt über eine Frage zum Stapelaustausch hinausgeht.

Barrycarter
quelle