Boundary-conforming invertible maps from the unit square onto the physical domain generate unfolded grids. On the other hand, if a map is not invertible, it produces folded grids and clearly it cannot be used in computations. The formulation of a variational principle which allows to separate the set of invertible maps from those maps that are not invertible is given.The functional, called the energy density of the map, depends on first derivatives of the coordinates and variable coefficients which form a symmetric positive-definite matrix. The variational principle states that the map is invertible if and only if it is the minimum of the energy density functional. It follows then, that the Euler equations of the functional describe all possible invertible (smooth and one-to-one) maps.Based on the energy density of grid deformation, a discrete analog of the variational principle is formulated and a complete proof is given. The energy density of grid deformation is a discrete functional with an infinite barrier on the boundary of the set of unfolded grids. The barrier property is very important for problems with moving boundaries and for moving adaptive grid methods because it ensures the generation of unfolded grids at every time step. Direct control of cell shapes, with prescribed cell width near the boundary, is used to achieve grid orthogonality on the boundary. Adaptation to the solution of the host equations is realized as in the adaptive-harmonic grid generator.The variational principle is implemented for the numerical solution of shallow water equations. Simulation of wind-induced currents in the Mojaiskii reservoir has been performed and results are presented. The adaptive numerical solution with boundary orthogonality increases the accuracy of the computations and is consistent with experimental data
