Problems "orbit"

Integrate (that is, solve) 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).
Consider the equation (in certain units) of equatorial motion of a planet around a star in General Relativity,
u(φ)'' + u(φ) = 1 + εu(φ)² ,
where u(φ) ≡ ¹/_r(φ) , r is the (reduced-circumference) 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,
```
y₀' = y₁
y₁' = 1-y₀+ε*y₀*y₀ 
```
- 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 
```
(Optional) Phase shifts.
The Schrödinger equation for the s-wave motion of a particle with mass m= in a short-range potential V(r) with energy E=ℏ²k²/(2m) is given as
```
-u''(r) + (2m/ℏ²)V(r)u(r) = k² u(r) ,
```
where u(r) is the radial wave-function.
For positive energies the boundary condition at r=0 is u(0)=0, u'(0)=1; and at r→∞ u(r)=Asin(kr+δ(k)), where A is a constant and δ(k) is the phase-shift.
Assume m=940MeV, V(r)=-50exp(-(r/1.7fm)²)MeV, ℏc=197.3 MeV*fm, and calculate the phase-shifts for the energies up to 100MeV by integrating the Schrödinger equation from r=0 to r≈10 and extracting the phase-shifts by matching the integrated solution to the asymptotic form.