Seepage with free surface is essentially one type of nonlinear free boundary problems. In order to solve this problem in global domain

the nodal virtual flux method deducts the virtual flux gradually in each iteration. It has the advantages of weak mesh dependency

fast convergence rate

etc. To reveal its internal theoretical foundation

this paper builds an equivalent bridge between the nodal virtual flux method and variational inequality of Signorini type through a complementarity constraint. The criterion of free surface is optimized by introducing an amplification coefficient for transition region

verified through a sand flume model test. The comparisons show that the improved algorithm has better numerical stability and higher curvature tolerance

and provides an effective approach for optimization design of seepage control structure with large-scale grid.