Planetensystem mit Mathematica

Dieses Notebook zeigt eine weitere Möglichkeit, einen Animations-Effekt zu erzielen.

Berechnung durchführen und Animation vorbereiten 
In der ersten Zelle werden die Newtonschen Differentialgleichungen für ein gravitatives System aus drei Objekten gelöst, von denen eines (rot dargestellt) eine vernachlässigbar kleine Masse besitzt ("restringiertes Dreikörperproblem"). Dabei bedeuten beispielsweise: m1 = Masse der Sonne, m2 = Masse eines schweren Planeten (etwa Jupiter), das dritte Objekt ist ein leichter Planet (etwa die Erde) oder ein Asteroid. Die Variablen xc und yc bezeichnen die Koordinaten des Schwerpunkts. Schließlich wird die Animation vorbereitet.

Clear[m1, m2, x1, y1, x2, y2, x3, y3, xc, 
        yc];
\(m1 = 0.8; \ m2 = 0.2; rad1 = m1/8; \ 
    rad2 = m2/8;
 (*\ \(Berechnung\)\(:\)\ *) \)
mred = (m1 m2)/(m1 + m2);
r[t_] = 
        Sqrt[(x1[t] - x2[t])^2 + (y1[t] - y2[t])^2];
r13[t_] = 
        Sqrt[(x1[t] - x3[t])^2 + (y1[t] - y3[t])^2];
r23[t_] = 
        Sqrt[(x2[t] - x3[t])^2 + (y2[t] - y3[t])^2];
\(loesungDgl = \(NDSolve[{\(x1''\)[t] == 
              m2*\((x2[t] - x1[t])\)/r[t]^3, \(y1''\)[t] == 
              m2*\((y2[t] - y1[t])\)/
                  r[t]^3, 
				\(x2''\)[t] == 
              m1*\((x1[t] - x2[t])\)/r[t]^3, \(y2''\)[t] == 
              m1*\((y1[t] - y2[t])\)/
                  r[t]^3, 
				\(x3''\)[t] == 
              m1*\((x1[t] - x3[t])\)/r13[t]^3 + 
                m2*\((x2[t] - x3[t])\)/
                    r23[t]^3, 
				\(y3''\)[t] == 
              m1*\((y1[t] - y3[t])\)/r13[t]^3 + 
                m2*\((y2[t] - y3[t])\)/
                    r23[t]^3, 
				x1[
                0] == \(-mred\), 
            y1[0] == 0, \(x1'\)[0] == 
              0, \(y1'\)[0] == \(-0.3\), 
				x2[
                0] == 1 - mred, 
            y2[0] == 0, \(x2'\)[0] == 0, \(y2'\)[0] == 
              0.71, 
\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ x3[0] == 
              0.5, \(x3'\)[0] == \(-0.71\)*Cos[Pi/6], 
            y3[0] == Sqrt[3]/2, \(y3'\)[0] == 
              0.71*Sin[Pi/6]}, 
				{x1[t], y1[t], 
            x2[t], y2[t], x3[t], y3[t]}, {t, 0, 12}]\)[\([1]\)]\)
x1[t_] = x1[t] /. loesungDgl;
x2[t_] = x2[t] /. loesungDgl;
x3[t_] = x3[t] /. loesungDgl;
y1[t_] = y1[t] /. loesungDgl;
y2[t_] = y2[t] /. loesungDgl;
y3[t_] = y3[t] /. loesungDgl;
xc[t_] = (x1[t] m1 + x2[t] m2)/(m1 + m2);
yc[
          t_] = (y1[t] m1 + y2[t] m2)/(m1 + 
              m2);
  
animieren := 
        Table[ParametricPlot[{{rad1 Cos[u] + x1[t] - xc[t], 
                rad1 Sin[u] + y1[t] - yc[t]}, {rad2 Cos[u] + x2[t] - xc[t], 
                rad2  Sin[u] + y2[t] - yc[t]}, {0.01 Cos[u] + x3[t] - 
                  xc[t], 0.01 Sin[u] + y3[t] - 
                  yc[t]}}
, {u, 0, 2*Pi}, 
            PlotRange -> {{-2, 2}, {-2, 2}}, 
            AspectRatio -> Automatic, 
            PlotStyle  ->  {RGBColor[0,  0, 0], RGBColor[0,  0, 0], 
                RGBColor[1,  0, 0]}], {t, 0, 12, 0.1}];

Animation ausführen 
Wurde die obige Zelle ausgeführt, so ist jetzt alles Nötige definiert, und Sie können die Animation starten:

animieren;