My code: Code (Python): import numpy as np import matplotlib.pylab as plt import math import random from scipy import integrate R1 = .001 R2 = 7 def G(r,theta): sigma = random.randint(4000., 7000.)/1000. # width of beam is 4 - 7mm r0 = random.randint(0, R1*1000.)/1000. #random centroid theta0 = random.uniform(0, np.pi*2) #angle return (100/(np.pi*2*sigma**2))*(np.e**(-((r**2 + r0**2 - 2*r*r0*(np.cos(theta)*\ np.cos(theta0) + np.sin(theta)*np.sin(theta0)))/2*sigma**2)))*r # this r is here because the integrand is dr and dtheta, NOT dx and dy any longer!!! #amplitude is 100 RTI = integrate.nquad(G, [[0,R2],[0,2*np.pi]]) #Integration print RTI I keep running this code, and everything seems to be working, but I keep getting results like: (0.1406360958428128, 0.012357271987783056) Also, the code itself takes quite a while to run and I'm sure this is an indication of something...just not sure where it's going on. As you can see, the error is almost the value of integration itself, and this error is to large. I've been looking at the code but am not sure why such a large error is arising. If anyone has any ideas, it would be greatly appreciated!