I am trying to solve an integral given as

$int_{-pi}^{pi} int_{-1}^{1} d(cos_theta) d(phi) sin_{theta}*cos_{phi}*g $

using scipy integrate.trapezoid method.

The code is as follows :

`def dist(a, vz):
     return np.exp(-(vz-1)**2)/(a)**2

 a = random.random()

 nvz = 256                               
 nph = 256                               

 vz = -1.0 + (0.5 + np.arange(nvz)) * (2.0 / nvz)
 ph = -np.pi + (0.5 + np.arange(nph)) * (2.0 * np.pi / nph)

 for i in range(nph):
     for j in range(nvz):
         g1 = dist(sigma_nu, vz[j])
         g.append(g1)

 g = np.array(g)
 dGamma = (2.0/nvz) * (2*np.pi / nph)
 g = g.reshape(nvz, nph)
 nnu = g.sum()*dGamma
 g = g/nnu

 dphi = (2*np.pi / nph)
 dtheta = (2.0/nvz)

 cos_theta = vz[np.newaxis,:]
 sin_theta = np.sqrt(1-cos_theta**2)
 phi = ph[:,np.newaxis]
 cos_phi = np.cos(phi)
 sin_phi = np.sin(phi)

 I = integrate.trapezoid((integrate.trapezoid(sin_theta*g*cos_phi, dx = dtheta, axis = 1)),  

 dx = dphi, axis = 0`

Now since ‘g’ is only dependent on theta, the integral ‘I’ should give 0 as the result. But it is not happening so. I am not able to figure out what is the mistake I am making.

Please suggest me a way out.