Binary dielectric metasurfaces are arrays of sub-wavelength structures that act as a thin layer of artificial material. They are generally lossless and relatively simple to fabricate since only a single structuring step is required. By carefully designing the metasurface, the phase, amplitude and the polarization of the incident light can be controlled at will. In practice, fabrication constraints and the limited choice of materials reduce what can be done with metasurfaces. But a wide range of functionalities can be implemented with the proper design techniques and knowledge. This thesis contributes to both. The modes are key to understand the phenomena occurring inside a metasurface. To facilitate the analysis of the modes, the Poynting operation, which is related to the Poynting vector, is introduced. We describe how this operation can be used to reformulate the boundary condition in order to estimate the reflection and transmission coefficients with reduced knowledge on the modes involved, and to orthonormalized the modes. The Fourier modal method, which is the method used for the rigorous simulation of metasurfaces in this work, is improved in order to facilitate the access to valuable information that can be used to better understand the phenomena occurring inside a metasurface, and to speed up the design and optimization process. This method computes the eigen-modes present in the metasurface. To better analyze them, the eigen-modes are orthonormalized using the Poynting operation. We show that most of the modes can be filter out in order to simulate a metasurface with different thicknesses in milliseconds. From the analysis of the modes propagating in metasurfaces, two types of metasurface are identified: single-mode metasurfaces and multi-mode metasurfaces. For single-mode metasurfaces, we provide design techniques that translates the desired response into internal properties of the metasurfaces. For multi-mode metasurfaces, the concept of self-coupling mode is developed. We show that, based on this concept, the angular and spectral response of a metasurface can be interpolated safely with a few simulations even if high-Q resonances are present. For both types of metasurfaces, examples of design are provided. Gradient-based optimization methods allow to obtain the optimal metasurface in a few iterations, but the derivative of the merit function is necessary. The adjoint method computes the functional derivative of the merit function with respect to the permettivity and permeability. We provide the equations of the adjoint method and apply them to diffractive optical elements such that they can be used in conjunction with the Fourier modal method. This thesis contributes to design challenges for complex electromagnetic problems,and it does not only provides concepts, but also the tools to put them into operation.