Exercise "orbit"

Objective

Learn to pass functions into subroutines; learn to use ode.rk23 to integrate ordinary differential equations (ODE).

Exercises

1. Integrate (that is, solve numerically) the ordinary differential equation,

   y'(x) = y(x)*(1-y(x))
   

   from x=0 to x=3 with the initial condition y(0)=0.5. Compare with the analytic result (should be the logistic function).

2. Consider the equation of equatorial motion (in certain units) of a planet around a star in General Relativity,

   u''(φ) + u(φ) = 1 + εu(φ)² ,

   where u(φ) ≡ ¹/_r(φ) , r is the (circumference-reduced) radial coordinate, φ is the azimuthal angle, ε is the relativistic correction (on the order of the star's Schwarzschild radius divided by the radius of the planet's orbit), and primes denote the derivative with respect to φ.

   1. Integrate this equation with ε=0 and initial conditions u(0)=1, u'(0)=0 : this should give a Newtonian circular motion.

   2. Integrate this equation with ε=0 and initial conditions u(0)=1, u'(0)≈-0.5 : this should give a Newtonian elliptical motion. Hint: u'(0) shouldn't bee too large or you will lose your planet.

   3. Integrate this equation with ε≈0.01 and initial conditions u(0)=1, u'(0)≈-0.5 : this should illustrate the relativistic precession of a planetary orbit.

   Hints:

   • You can rewrite this second order equation as a system of ordinary differential equations (ODE) by introducing the functions y₀=u and y₁=u', which gives the following ODE,

     y0' = y1
     y1' = 1-y0+ε*y0*y0 

   • If you have calculated the data-file as two columns with φ and u, you can plot the orbit with the following gnuplot/pyxplot command,

     plot "data" using (1/$2)*sin($1):(1/$2)*cos($1) with lines notitle