HANDS ON for COMPUTATIONAL METHODS IN COSMOLOGY

• Hands on Preparation:
  units, variables, classes, loops, first program and plotting)
• Hands on I (90 min):
  Numerical Integration (Euler and Leapfrog)
• Hands on II (90 min):
  Treating N-Bodies (Central potential and self-gravity)

Hands On preparation

First step: units and acceleration; a first program

• How do you re-define natural constants when you choose a different unit system?
  Example: the gravitational constant is 6.674 × 10⁻⁸ cm³ g⁻¹ s⁻² in the cgs system and 6.674 × 10⁻¹¹ m³ kg⁻¹ s⁻² in the SI system.
• In your first program, create variables to define a unit system and calculate the value for the gravitational constant G in these units.
• Calculate the different contributions of the gravitational acceleration at the Earth’s surface by the Earth itself, the Moon, the Sun, the Galaxy, and the Virgo cluster. What is the gravitational acceleration on the surface of a neutron star? What is the acceleration Pluto feels from the Sun? What is the acceleration the Voyager spacecraft feels from the Sun? What is the acceleration a galaxy feels at the border of a massive galaxy cluster? Check how accurately Newton’s law of gravity can currently be measured on Earth.

Second step: gravitational forces; vectors and loops; gnuplot

• Check what is the Lense-Thirring precession and how it compare to Newtonian acceleration?
  Hint: have a look at equations (97,98) in this paper.
• Write a program that outputs a comparison between the Newtonian acceleration and the Lense-Thirring precession of a star close to the galactic centre. Create a loop and calculate these values for a series of distances and plot the output with gnuplot.

Hands on I: Numerical integration (Euler and Leapfrog

1. • Write a program that just integrates the earth orbit around the sun using a forward Euler scheme: df = f' dt
     It should output the time and the x,y,z coordinates for each timestep.
     Can you write it using 'natural' units?
     
   • Make the program compute the total energy and print this out as well.
     How does the system evolve in time for different timestep sizes?
     How does it evolve if you change the orbit to be more and more eccentric?
     
     Try to define a structure which collects all the variables you need for describing the system at a given time.
     Try to define sub-routines where you do the individual steps.
     Try to pass a reference to your structure to these sub-routines.
2. • Modify your program by implementing a more advanced integrator.
     Implement both leap frog variants and give a possibility to choose the integrator.
     Repeat the same experiments as in 1) (varying timestep and eccentricity).
     Modify the program to output only every nth orbit (or change your plotting routine to plot only every nth orbit) in order to better see the evolution of the orbits.
     
   • Allow the timestep to be variable.
     Can you construct a useful criterion for setting the timestep size?
     What does the eccentric orbit look like now?

Hands on II: Treating N-Body systems

1. • Take your leap-frog integrator from last time and extend the program to handle multiple particles.
     Setup a solar system and plot the orbits.
     How do you implement your variable timestep criteria in this case?
     Setup the whole solar system as test case.
     Plot the orbits of the planets in gnuplot!
     
   • Change the central potential to be a plummer potenial.
     Include as output the effective potential for each particle. How does it look like?
     Distribute randomly particles (maybe some tens ?) with zero initial velocity and let the system evolve.
     Watch the evolution of the particles in phase space
     
2. • If we have time: Instead of a potential implement self gravity between all particles.
     For example add the sun and rerun the solar system.
     
     Try to read the initial conditions directly from a file.
     Can you give the file name from the command line top the program ?