The Fermi Gamma-ray Space Telescope greatly increased the number of known pulsars in the gamma-ray band. It was discovered that pulsars often emit a significant fraction of their spin-down energy in gamma-rays, offering a new powerful diagnostic tool of neutron star magnetic field structure. The modeling of gamma-ray data suggests that in order to interpret observations one needs to understand the field geometry and the plasma state in the emission region. In recent years, significant progress has been achieved in understanding the magnetospheric structure in the limit of abundant plasma supply. However, the very presence of dense plasma everywhere in the magnetosphere is not obvious. The region where the observed emission is produced still remains a mystery. To address this from first principles, a full kinetic treatment, which includes the physics of plasma production and particle acceleration, is necessary. In this proposal we will perform kinetic simulations of magnetospheres using relativistic particle-in-cell codes. We will construct both gamma-ray light curves and spectra with minimal assumptions. This will allow meaningful interpretation of Fermi observations of pulsars and answer long-standing puzzles about the physics of pulsar mechanism. These investigations are consistent with NASA's strategic goal of discovering the origin, structure, and evolution of the universe.