<p> In this paper, we propose a novel near-field beamforming (BF) design for a Large Intelligent Surface (LIS) implemented as a discretized two-dimensional (2D) antenna array. We firstly investigate the definitions of near-field and far-field regions of the LIS, and determine the Fraunhofer distance which scales up linearly in the surface-area. Then, we derive the Fresnel near-field region where variations of amplitudes and angles are both negligible, as long as the distance from a user-equipment (UE) to the LIS is larger than a threshold that scales up linearly in the diameter of the LIS. Hence, in the majority region of near-field, only phase variations worsen the quality of received signal and yield significant losses of array gains. Following this observation, we present a decomposition theorem proving that in the Fresnel region, the optimal three-dimensional (3D) near-field BF can be decomposed into a conventional 2D far-field BF and a one-dimensional (1D) near-field BF. The 2D far-field BF compensates phase variations caused by the azimuth and elevation angles, while the 1D near-field BF compensates phases variations caused by distance differences among the UE and antenna-elements on the LIS. Such a 2D+1D BF design is fully compatible with existing 2D far-field BF schemes, and the 1D near-field BF can be treated as an additional step. From this perspective, we further derive an optimal codebook design of the 1D near-field BF with Lloyd-Max method. Numerical results show that the proposed 2D+1D BF is effective to recover losses of array gains in the near-field, and also robust against distance mismatches. Moreover, with a few optimized quantization points in the codebook, the 1D near-field BF performs close to optimal. </p>