A variational ground state of the repulsive Hubbard model on a square lattice is investigated numerically for an intermediate coupling strength (U = 8t) and for moderate sizes (from 6 × 6 to 10 × 10). Our ansatz is clearly superior to other widely used variational wave functions. The results for order parameters and correlation functions provide new insight for the antiferromagnetic state at half filling as well as strong evidence for a superconducting phase away from half filling.PACS numbers: 71.10. Fd,74.20.Mn, The Hubbard model plays a key role in the analysis of correlated electron systems, and it is widely used for describing quantum antiferromagnetism, the Mott metalinsulator transition and, ever since Anderson's suggestion [1], superconductivity in the layered cuprates. Several approximate techniques have been developed to determine the various phases of the two-dimensional Hubbard model. For very weak coupling, the perturbative Renormalization Group extracts the dominant instabilities in an unbiased way, namely antiferromagnetism at half filling and d-wave superconductivity for moderate doping [2,3]. Quantum Monte Carlo simulations have been successful in extracting the antiferromagnetic correlations at half filling [4,5], but in the presence of holes the numerical procedure is plagued by the fermionic minus sign problem [6]. This problem appears to be less severe in dynamic cluster Monte Carlo simulations, which exhibit a clear tendency towards d-wave superconductivity for intermediate values of U [7].Variational techniques address directly the ground state and thus offer an alternative to quantum Monte Carlo simulations, which are limited to relatively high temperatures. Previous variational wave functions include mean-field trial states from which configurations with doubly occupied sites are either completely eliminated (full Gutzwiller projection) [8,9,10] or at least partially suppressed [11]. Recently, more sophisticated wave functions have been proposed, which include, besides the Gutzwiller projector, non-local operators related to charge and spin densities [12,13]. Our own variational wave function is based on the idea that for intermediate values of U the best ground state is a compromise between the conflicting requirements of low potential energy (small double occupancy) and low kinetic energy (delocalization). It is known that the addition of an operator involving the kinetic energy yields an order of magnitude improvement of the ground state energy with respect to a wave function with a Gutzwiller projector alone [14]. In this letter, we show that such an additional term allows us to draw an appealing picture of the ground state, both at half filling and as a function of doping (some preliminary results have been published [15,16]).In its most simple form, the 2D Hubbard model is composed of two terms,Ĥ = tT + UD , witĥHere c † iσ creates an electron at site i with spin σ, the summation is restricted to nearest-neighbor sites and n iσ = c † iσ c iσ . We consider a square lattice with peri...