You can get rid of the warning with the following replacement:
#ds=Measure('ds')[boundary_parts]
ds=Measure('ds',subdomain_data=boundary_parts)
The fact that the flux at the bottom exactly matches the flux at the top is not surprising, since the exact solution just varies linearly from the top of the domain to the bottom. This solution is in your discrete space, so you will recover it exactly with Galerkin’s method. For some additional information on computing conservative fluxes on Dirichlet boundaries in more complicated scenarios, see the discussions here and here.