ABSOLUTE AND RELATIVE MOTION OF A BINARY STAR SYSTEM
K 7.3

Bernhard Schnizer

Institut fuer Theoretische Physik
Technische Universitaet Graz
A-8010  Graz, Austria


2002-05-03

Start Button for Animations

start = 
        ButtonBox[
          StyleBox["Run Animation", "Title", 
            FontFamily -> "Helvetica"], 
          ButtonFunction :> 
            FrontEndExecute[{SelectionMove[EvaluationNotebook[], Before, 
                    ButtonCell] SelectionMove[EvaluationNotebook[], 
                    Previous, CellGroup] SelectionAnimate[
                    SelectedNotebook[], AnimationDisplayTime -> zeit]}], 
          ButtonEvaluator -> Automatic, 
          Active -> True];

Equations of motion for relative distance

Remove["@"]

Equations of motion for the relative coordinate  r :

r  =  {x[t], y[t]};
v  =  D[r, t];
b  =  D[v, t];

rq  =  Sqrt[r . r]^3

Total mass is  m1 + m2 = 1; gravitational constant = 1. Then the equation of motion
for the relative distance  r  is :

sy0  =  Thread[b  ==  - r/rq]

The initial data for the two masses give those,  r0 and v0 , for relative distance  r.
These two 2-dim. are extended to 3-dim. vectors,  rr  und  vv, so that the 
vector products may be computed:

r0  =  r20  -  r10
v0  =  v20  -  v10
rr  =  Join[r0, {0}];
vv  =  Join[v0, {0}];

Angular momentum for the relative coordinate is:

lo  =  Cross[rr, vv]

ll  =  lo[[3]]

Total energy of the relative coordinate is :

e0  =  v0 . v0/2  -  1/Sqrt[r0 . r0]

The Laplace Lenz vector for the relative coordinate is:

ao  =  rr/Sqrt[rr . rr]  -  Cross[vv, lo]

The equation of the trajectory of the relative coordinate is:

bg[ph_]  =  ll^2/(1  +  ao[[1]] Cos[ph])

pg  =  
        ParametricPlot[-bg[ph] {Cos[ph], Sin[ph]}, {ph, 0, 2  Pi},  
          AspectRatio  ->  Automatic];

The distance of the pericenters is:

bg[0]

The distance of the apocenters is:

bg[Pi]

a2  =  %  +  %%

Large semi-axis :

aa  = a2 /2

Excentricity :

e  =  Abs[ aa ao[[1]]]

Small semi-axis :

bb  =  Sqrt[aa^2  -  e^2]

c  =  bg[Pi] -  aa;  
    vl  =  Line[{{c, bb}, {c, -bb}}];
bp  =  {PointSize[ .02], Point[{0, 0}], 
          Point[{c + e, 0}]};

The elliptic trajectory of the relative coordinate together with the focal points:

prr  =  Show[pg, Graphics[vl], Epilog  ->  bp];

Uniform motion of the center of gravity:

rso  =  m1 r10  +  m2 r20
vso  =  m1 v10  +  m2 v20

rs[t_]  =  rso  +  vso t

Integrating the equations for relative motion

Motion of the center of gravity is parametrized with the angle ph; the true anomaly with time  t.
Die Schwerpunktsbewegung ist mit der Zeit  t  parametrisiert, die obige Bahngleichung bg[ph]
These two are connected by the angular momentum:

					L   =   m (r^2) d(ϕ/dt)

 We do not integrate this equation. In place of it the equation of motion for the relative distance  r .

an0  =   {Thread[(r  /.  t  ->  0)  ==  
            r0],   
             Thread[(v  /.  t  ->  0)  ==  v0]}

dg0  =  Join[sy0, an0] // Flatten

tmax  =  30;

s0  =  NDSolve[dg0, {x, y}, {t, 0, tmax}];

A plot of this solution again gives the trajectory of relative distance:

p0  =  
        ParametricPlot[Evaluate[{x[t], y[t]} /. s0], {t, 0, tmax},  
          AspectRatio  ->  Automatic];

Show[prr];

Total motion of the binaries

The motion of the center of gravity augmented by the relative motion weighted with the 
masses m1, m2 respectively gives the motion of the masses   m1   und   m2 .

ct  =  AbsoluteThickness[2];

p1  =   ParametricPlot[
          Evaluate[(rs[t]  -  m2 {x[t], y[t]}) /. s0], {t, 0, 
            tmax}, 
PlotRange  ->  All,  
          PlotStyle  ->  {{ct, RGBColor[ .9, 0, 0]}}];

p2  =   ParametricPlot[
          Evaluate[(rs[t]  +  m1 {x[t], y[t]}) /. s0], {t, 0, 
            tmax}, 
PlotRange  ->  All,  
          PlotStyle  ->  {{ct, RGBColor[0, 1, 0]}}];

Motion of the center of gravity is:

ps  =  
        ParametricPlot[
          Evaluate[rs[t]], {t, 0, 
            tmax}, 
PlotStyle  ->  {{ct,  
                GrayLevel[ .7]}}];

So total motion of the binaries is:

tp  =  10;  sp  =  Point[rs[tp]]; 
po1  =  
        Point [Flatten[(rs[tp]  -  m2 {x[tp], y[tp]}) /. 
              s0]];
po2  =  
        Point [Flatten[(rs[tp]  +  m1 {x[tp], y[tp]}) /. 
              s0]];

fp  =  
        Show[p1, p2, ps,  
          Epilog  ->  {PointSize[ .01], sp, po1, po2}];

For a  center of gravity at rest  the relative motion of the binaries is:

\(\(pr1\  = \ 
        ParametricPlot[
          Evaluate[\(-\ m2\)\ {x[t], y[t]}\  /. \ s0], 
{t, 0, 
            tmax}, \ \ \ AspectRatio\  -> \ Automatic, \ 
          PlotStyle\  -> \ {{ct, RGBColor[ .9, 0, 0]}}];\)\)

\(\(pr2\  = \ 
        ParametricPlot[
          Evaluate[m1\ {x[t], y[t]}\  /. \ s0], 
{t, 0, 
            tmax}, \ \ \ AspectRatio\  -> \ Automatic, \ 
          PlotStyle\  -> \ {{ct, RGBColor[0, 1, 0]}}];\)\)

pq0  =   Point[{0, 0}];
pq1  =  
        Point [Flatten[(- m2 {x[tp], y[tp]}) /. 
              s0]];
pq2  =  
        Point [Flatten[(   m1 {x[tp], y[tp]}) /. s0]];

fq  =  
        Show[pr1, pr2,  
          Epilog  -> {PointSize[ .01], pq0,  pq1, pq2},  
          AspectRatio  ->  Automatic];

Print["         m1 = ", m1, "         m2 = ", 
      m2]
Show[fp];  Show[fq];

Relative motion, Time points

pq0  =   Point[{0, 0}];

np  =  12;  
    tk  =  Sqrt[4 π^2 aa^3];
lp  =  Table[tp, {tp, 0, tk,  tk/np}];

lrp  =  
        Table[Evaluate[Flatten[{x[lp[[k]]], y[lp[[k]]]} /. s0]] // 
            Chop, {k, Length[lp]}];

prp1  =  
        Graphics[Map[Point, - m2 lrp]];
prp2  =  Graphics[Map[Point,    m1 lrp]];

Show[pr1, pr2, prp1, prp2,  Epilog  ->  {pq0},  
        ImageSize  ->  500];

Simulation relative motion

pq0  =   Point[{0, 0}];

np  =  40;  
    tk  =  Sqrt[4 π^2 aa^3];
lp  =  
        Table[tp, {tp, 0, tk,  tk/np}];
lp2   =  
        Table[tp, {tp, tk, 2  tk,  tk/np}];
lp3  =  Table[tp, {tp, 2  tk, 3  tk,  tk/np}];

lrp  =  
        Table[Evaluate[Flatten[{x[lp[[k]]], y[lp[[k]]]} /. s0]] // 
            Chop, {k, Length[lp]}];

prp1  =  Map[Point, - m2 lrp];
prp2  =  Map[Point,    m1 lrp];

sr = Table[Show[pr1,pr2,Graphics[prp1[[k]]],Graphics[prp2[[k]]], Epilog -> {pq0}, ImageSize -> 500],{k,Length[lrp]}];

sp  = rs[lp] // Transpose ;  sp2  = rs[lp2] // Transpose ;  
    sp3  = rs[lp3] // Transpose ;
lsp  =  Map[Point, Join[sp, sp2, sp3]];

sm  =  Show[Graphics[Join[{GrayLevel[ .9]}, lsp]]];

l11  =  sp  -  m2 lrp;  l12  =  sp2  -  m2 lrp;  
    l13  =  sp3  -  m2 lrp;
l21  =  sp  +  m1 lrp;  l22  =  sp2  +  m1 lrp; 
    l23  =  sp3  +  m1 lrp;

pqp1  =  
        Map[Point, Join[l11, l12, l13]];
pqp2  =  Map[Point, Join[l21, l22, l23]];

Show[Graphics[pqp1], Graphics[pqp2], 
        ImageSize -> 500];

Show[sm, p1, p2, Graphics[pqp1], Graphics[pqp2]];