These strains are a lot smaller than \epsilon_{xx}, and the error decreases with refinement. In general, you expect some discretization error when your displacement function space does not contain the exact solution. In the case of pure bending, your exact solution involves quadratic variation of the vertical displacement u_y. You can see that the error in \epsilon_{xy} and \epsilon_{yy} becomes negligible when you use quadratic elements for the displacement:
#V = VectorElement("Lagrange", mesh.ufl_cell(), 1)
V = VectorElement("Lagrange", mesh.ufl_cell(), 2)