arXiv:2611.08276v1 Announce Type: new Abstract: A hybrid virtual element formulation (HVEM) is presented for transient flow in variably saturated soils governed by Richards' equation. The nodal heads are the only global degrees of freedom, and a Darcy flux is assumed in each polygonal element and condensed. The flux combines a divergence-free part with a particular field of linear divergence that absorbs storage and sources, so that the water balance of every element holds identically, at every iterate and for any time step, and that of the whole domain at convergence. To guarantee rank-sufficiency with no stabilisation parameter, the divergence-free fields are harmonic polynomials ordered by their behaviour under rotation, giving a closed-form rule for the number of flux parameters, derived on regular polygons and verified numerically on distorted ones and on elements with hanging nodes. Darcy's law is written in each element for the Kirchhoff potential, which removes the conductivity from the diffusion while the head stays continuous across material interfaces. A Kirchhoff mean of the conductivity makes hydrostatic equilibrium and gravity drainage exact on any mesh. The storage is lumped on the nodes with weights supplied by the assumed flux. All the element operators reduce to boundary integrals, so that the nonlinear iterations need only nodal values. The formulation is assessed on solutions in closed form, on published benchmarks and on rainfall over homogeneous and layered slopes, where independent solutions provide a further check. It converges at second order in the head on polygonal and locally refined meshes, also where gravity dominates, and closes the water budget in every test. On a layered slope, cut grids refined in a band along the surface reach the accuracy of a uniform grid twice as fine with about 40% of the nodes. Its structure, shared with the HVEM for elastoplasticity, suits coupled analyses of slopes.