This paper presents a weak Galerkin finite element method based on H(div) virtual element for Darcy flow on polytopal meshes. Specifically, the pressure is approximated by constant functions on both the polytopal elements and the element sides, respectively. Their discrete weak gradients are calculated in the local H(div) virtual element spaces whose basis functions never need to be explicitly computed. The method works on very general polytopal meshes. After obtaining the numerical pressure, the numerical velocity that is locally mass-conservative and has continuous normal fluxes is easily obtained by a simple L2 projection. The error estimates for the numerical pressure and numerical velocity are proved. The assembly of element stiffness matrix is introduced. Numerical experiments in two and three dimensions are reported to illustrate the good performance of the method.