Motivated by recent realizations of Dy2Ti2O7 and Ho2Ti2O7 spin ice thin films, and more generally by the physics of confined gauge fields, we study a model of spin ice thin film with surfaces perpendicular to the [001] cubic axis. The resulting open boundaries make half of the bonds on the interfaces inequivalent. By tuning the strength of these inequivalent "orphan" bonds, dipolar interactions induce a surface ordering equivalent to a two-dimensional crystallization of magnetic surface charges. This surface ordering can also be expected on the surfaces of bulk crystals. In analogy with partial wetting in soft matter, spins just below the surface are more correlated than in the bulk, but not ordered. For ultrathin films made of one cubic unit cell, once the surfaces are ordered, a square ice phase is stabilized over a finite temperature window, as confirmed by its entropy and the presence of pinch points in the structure factor. Ultimately, the square ice degeneracy is lifted at lower temperature and the system orders in analogy with the well-known F -transition of the 6-vertex model.Highly frustrated magnets have been shown to host an astonishing array of exotic many-body phenomena, taking us far from the conventional paradigms of collective magnetic behavior [1]. The formulation of the local frustrated constraints in terms of effective gauge fields has revolutionized our perspective on these systems in both the classical and quantum domains. Depending on the system, this gauge symmetry can take the form of electromagnetism [2-5] with photon and magnetic-monopole excitations [5][6][7][8][9][10][11][12][13] [28][29][30][31][32][33][34].In this paper, we study the dipolar spin ice model in thin film (slab) geometry (slab-DSI) defined by Eq. (1) below. As illustrated in Fig. 1, the slab geometry renders the nearest-neighbor bonds on the surface inequivalent: either a bond belongs to a "bulk tetrahedron", or it is an "orphan bond", belonging to a "virtual tetrahedron" which has been "cleaved away" at the interface. Each virtual tetrahedron may therefore carry a magnetic surface charge that can propagate to and from the surface into the bulk. The orphan bonds act as effective chemical potentials for surface charges, allowing the slab-DSI to go through a surface charge-ordering transition. This is monopole crystallization [35] in two-dimensions, driven by the magnetic Coulomb potential between surface charges. Most notably, the crystallization is limited to the microscopic surface layer with no penetration into the bulk. The thinnest slab where this phenomenology can be explored systematically is one cubic unit-cell thick, containing three layers of tetrahedra. Below the surface ordering temperature at T so , the central layer emerges in the form of a constrained square ice system, which supports a state with extensive degeneracy described by a two-dimensional Coulomb phase [36]. Dipolar interactions ultimately break the symmetry of the Coulomb phase in the same way as in the venerable antiferroelectric F -model [37]...