Dear Nathan,

Suppose I have this function (given below), where ne, and Temperature is the value corresponding to the grid, but this will be a yt-field, number density, and temperature of the respective grid.

def max_jut_dist(p, ne, Temperature) :

func = ne.d / (4.0*np.pi * (m_electron.d * speed_of_light.d)**3)
theta = boltzmann.d * Temperature.d / (m_electron.d * speed_of_light.d**2)
gamma = np.sqrt(1.0 + (p / (m_electron.d * speed_of_light.d))**2)
func /= theta
func /= syspc.kn(2, 1.0/theta)
func *= np.exp( -gamma/theta)
return func

if I call this function inside the field to compute normalization of this function on each grid over p variable:

using ---->>> scipy.integrate.quad( max_jut_dist, g1, g2, args=( data["electron_density"] , data["Temperature"]) )

but this will give an error, shape () doesn't match the broadcast shape (16, 16, 16)

I hope this will explain what I want to do.

Thanks

-Prateek