Disagreement between Analytic Convolution and FFT Convolution
After playing with the notebook for a while, I figured out the following:


The result from scipy is shifted by one pixel in each axis (probably due to its definition of the "same" mode). Or, alternatively, both the analytic result and pyfftw results are shifted the opposite way.

After working around that, the result is exactly (within error) the same as the pyfftw result. I did this as follows (though this drops some data on the edges):




B_padded = np.pad(B, ((0, 1), (0, 1)), 'edge')[1:, 1:]
C1 = scipy.signal.convolve(A, B_padded, 'same', 'direct') * delta**2



The analytic solution (which peaks at 4.0) is not quite right for your input data (31x31 pixels at 1.0, edges at 0.5, and corners at 0.25). The peak value (where the image and the window align perfectly) should be:




box_bulk = 1.0**2 * (31 * 31)
box_edges = 0.5**2 * 31 * 4
box_corners = 0.25**2 * 4
peak = (box_bulk + box_edges + box_corners) / 1024 * 4  # ~= 3.876


which is the answer produced by both scipy and pyfftw.

Gist (also includes plotting tweaks).

(disclaimer: friend of asker)