*** Interpolating numerical data *** May 22, 2010: E.Komatsu In many cases, analytical functions are not available and one must interpolate numerical data. In such cases, it is useful to have an "interpolating function" which can be used as a function but all it does is to interpolate the numerical data. Here we provide a sample routine for interpolating the precomputped linear matter power spectrum at z=0. *** P(k) is in units of h^-3 Mpc^3, and k is in units of h Mpc^-1 *** (This is the same as CAMB definition.) For P(k), we provide the sample data, "wmap5baosn_max_likelihood_matterpower.dat," which was generated using CAMB code for the maximum likelihood parameters given in Table I of Komatsu et al.(2008) [WMAP 5-year interpretation paper] with "WMAP5+BAO+SN". The input file for CAMB is also provided (wmap5baosn_max_likelihood_params.ini). NOTE THAT THIS POWER SPECTRUM IS COMPUTED AT Z=0. (1) Put "USE linearpk" at the beginning of the program (2) Define the double-precision function "linear_pk". [double precision :: linear_pk] (3) Declare that it is an external function. [external linear_pk] (4) CALL open_linearpk(filename,n) where "filename" contains the name of power spectrum file (character), and n is the number of lines in the file (integer) (5) Then, call it as a usual function, "linear_pk(ak)" where "ak" is the wavenumber. The argument is in double precision. (6) CALL close_linearpk when you are done. - ak: k in units of h Mpc^-1 - linear_pk: P(k) in units of h^-3 Mpc^3 - To compile and use the sample programs (given below), edit Makefile and simply "./make" - It will generate executables called "sample" and "compute_linearpk" ========================================================================= A simple program to use "linearpk.f90" (also included as "sample.f90" in this directory): USE linearpk double precision :: k_ov_h,linear_pk character(len=128) :: filename integer :: n external linear_pk filename='wmap5baosn_max_likelihood_matterpower.dat' n=896 ! # of lines in the file CALL open_linearpk(filename,n) k_ov_h=0.01d0 ! h Mpc^-1 print*,'P(k) at k=',k_ov_h,' h Mpc^-1 is',linear_pk(k_ov_h),' h^-3 Mpc^3' CALL close_linearpk end A program for generating P(k) as a function of k is provided as "compute_linearpk.f90". The input linear spectrum is "wmap5baosn_max_likelihood_matterpower.dat". It looks like this: PROGRAM Compute_LinearPk USE linearpk ! A sample program for interpolating the numerical table of ! the linear matter power spectrum. ! - k is in units of h Mpc^-1 ! - P(k) is in units of h^-3 Mpc^3 ! May 22, 2010: E.Komatsu IMPLICIT none integer :: j double precision :: ak,linear_pk,dlnk,lnk character(len=128) :: filename integer :: n external linear_pk ! read in and tabulate P(k) filename='wmap5baosn_max_likelihood_matterpower.dat' n=896 ! # of lines in the file CALL open_linearpk(filename,n) ! now output P(k) as a function of k open(1,file='wavenumber_pk.txt') ak=4d-3 dlnk=1d-2 do while (ak<0.45d0) write(1,'(2E18.8)')ak,linear_pk(ak) lnk=dlog(ak)+dlnk ak=dexp(lnk) enddo dlnk=1d-1 do while (ak<15d0) write(1,'(2E18.8)')ak,linear_pk(ak) lnk=dlog(ak)+dlnk ak=dexp(lnk) enddo close(1) CALL close_linearpk END PROGRAM Compute_LinearPk