I am trying to find the expression for the rate of change of comoving volume with respect to redshift, that is $\frac{\mathrm{d}V_c}{\mathrm{d}z}$. In this paper (Hogg, David W.), the comoving volume differential is given as follows:
$$\mathrm{d}V_c=D_\mathrm{H}\frac{(1+z)^2D_\mathrm{A}^2}{E(z)}\mathrm{d}\Omega\mathrm{d}z$$
$$\implies\frac{\mathrm{d}V_c}{\mathrm{d}z}=D_\mathrm{H}\frac{(1+z)^2D_\mathrm{A}^2}{E(z)}\mathrm{d}\Omega$$
Where $D_\mathrm{H}$ is the Hubble distance, $D_\mathrm{A}$ is the angular diameter distance. These expressions can be expressed as follows:
$$D_\mathrm{H}=\frac{c}{H_0}$$
Where $c$ is the speed of light in vacuum, and $H_0$ is the Hubble's constant.
$$D_\mathrm{A}=\frac{D_\mathrm{M}}{1+z}$$
Where $D_\mathrm{M}$ is the transverse comoving distance given by:
$$D_\mathrm{M}=D_\mathrm{C},\textrm{ }\Omega_k=0$$
Where $\Omega_k=0$ implies spatial curvature 0, and $D_\mathrm{C}$ is the comoving distance given by:
$$D_\mathrm{C}=D_\mathrm{H}\int_0^z\frac{\,\mathrm{d}z'}{E(z')}$$
Now, this derivation can be substituted in $\frac{\mathrm{d}V_c}{\mathrm{d}z}=4\pi D_\mathrm{H}\frac{(1+z)^2D_\mathrm{A}^2}{E(z)}$, where $4\pi$ substitutes the solid angle element $\mathrm{d}\Omega$. Now, I want to plot this result as a function of redshift, for which I have written the following code in Python:
import numpy as np
import matplotlib.pyplot as plt
from astropy.cosmology import FlatLambdaCDM
cosmopar = FlatLambdaCDM(H0 = 67.8, Om0 = 0.3) #Defining the Cosmological model
speed_light = 3*10**5 #km/s
D_H = (speed_light/cosmopar.H(0).value)/10**3 #Gpc
Omega_m = 0.3
Omega_k = 0
Omega_Lambda = 1-Omega_m
def E(z_m): #Redshift z_m
return math.sqrt((1+z_m)^2*Omega_m+(1+z_m)^2*Omega_k+Omega_Lambda)
def Comov_vol(z_m): #in Gpc^3
return 4*math.pi*D_H*(((1+z_m)^2*(cosmopar.angular_diameter_distance(z_m).value/10**3)^2)/(E(z_m)))
#Plotting dVc_dz as a function of z:
ztrial = np.linspace(0,10,100) #Creating an array of redshift
comovol_arr = np.array([Comov_vol(v) for v in ztrial]) #Creating an array for corresponding comoving volume per redshift
plt.xlabel(r'Redshift')
plt.ylabel(r'$\frac{\mathrm{d}V_c}{\mathrm{d}z}$')
plt.plot(ztrial, comovol_arr)
This gives me the following plot:
Now, as one might notice, the magnitude of $\frac{\mathrm{d}V_c}{\mathrm{d}z}$ (in Gpc$^3$) quite large as compared to literature (Karathanasis, Christos, et. al.). I am referring to this paper for a project I am working on and the comoving volume per unit redshift vs redshift graph provided in this paper is as follows:
This graph is again given in Gpc$^3$, and the shape of both my plot and this literature plot is the same; however, the magnitude of $\frac{\mathrm{d}V_c}{\mathrm{d}z}$ seems to differ vastly. Can someone please help me understand where I am going wrong?
For context: I am new to both Astronomy and Coding. As an undergraduate Physics student, I am working on a project as a part of my summer internship and the paper which I refer to for the literature plot is the one required by my mentor, thank you.