Q:

Why does scipy provide a nonzero result for the integral where the result should be zero?

Given the following function:
d = y-a*sin(b*x-c)

I integrate this function over specific bounds for a, b and c using nquad from scipy:
def Iper(x,y,alb,aub,blb,bub,clb,cub):
    return nquad(lambda a,b,c: y-a*np.sin(b*x-c), [[alb,aub],[blb,bub], [clb,cub]], full_output=True)

using the boundaries below:
print(Iper(5.0,3.0,
           1.5,2.0,
           4.3,4.5,
           7.5,8.0))

I get the following result:
(0.09107123193614816, 1.2123921401960992e-14, {'neval': 9261})

However the result should be (best to my knowledge) zero. Plotting the resulting d values for the given bounds yields an area which is by no means intersected by the function of d. Below you see the function of d with the axis b, c and d. The green and violet curve are for the a values a = 1.5 and a = 2.0. The eight points represent the 'box'/area encompassed by the bounds of a, b, c and d... none of the two drawn functions intersect with this box area. How does it come that we have here an integral which does not equal zero? 

A:

I think you're misinterpreting what a 3 dimensional multivariate integral measures. It's not the 3d volume of the intersection of the domain and the curves (this is what a 2d multivariate integral measures), but the 4 dimensional hypervolume.

Looking at this integral step by step, first integrating with respect to a:

Then b:

Then c:

We see the integral is non-zero for the given values of x and y.

Note: when I evaluated this I got:
(0.06749972471060735, 1.031215181499445e-14, {'neval': 9261})

instead of your answer.

