*** Proper Angular Diameter Distance in units of h^-1 Mpc *** August 23, 2008: E.Komatsu Here we provide a simple subroutine, da.f90, for computing the proper angular diameter distance, Da(z), in units of h^-1 Mpc. IMPORTANT: the radiation density is ignored in the calculation, and therefore this routine cannot be used to compute the angular diameter distance to the decoupling epoch, z=1090, as radiation density should not be ignored there. For the other applications this approximation should be OK. (1) Put "USE cosmo" at the beginning of the program (2) Put "USE angular_distance" at the beginning of the program (3) Specify three cosmological parameters: - matter density parameter, om0 - dark energy density parameter, ode0 - dark energy equation of state parameter, w (4) CALL setup_da (5) You are now ready to use a double-precision function "da(z)" where z is also in double precision. - To compile and use the sample programs (given below), edit Makefile and simply "./make" - It will generate executables called "sample" and "compute_da" ========================================================================= A simple program to use "da.f90" (also included as "sample.f90" in this directory): USE cosmo USE angular_distance double precision :: da,redshift external da ! Specify three cosmological parameters in double precision. ! The data type has been defined in MODULE cosmo. ode0=0.723d0 om0=0.277d0 w = -1d0 CALL setup_da redshift=10d0 print*,'omega matter=',om0 print*,'omega de=',ode0 print*,'w=',w print*,'Da(z) at z=',redshift,'is',da(redshift),' h^-1 Mpc' end Another program, "compute_da.f90," for generating Da(z) over many redshifts is PROGRAM Compute_da USE cosmo USE angular_distance ! A sample program for computing the proper angular diameter distance ! in units of h^-1 Mpc ! August 23, 2008: E.Komatsu IMPLICIT none integer :: j double precision :: da,z external da ! specify four cosmological parameters ode0=0.723d0 om0=0.277d0 w=-1d0 ! tabulate da(z) by calling "setup_da" call setup_da ! now output da(z) as a function of z... open(1,file='redshift_da.txt') z=0d0 do j=1,10000 write(1,'(2F18.8)')z,da(z) z = z+0.01d0 enddo close(1) stop END PROGRAM Compute_da