In this article, both numerical and analytical studies of dielectric and orientational relaxation in a simple model system are presented. The system consists of point dipoles fixed at the sites of a simple cubic lattice. The dipoles are free to undergo rotational Brownian motion in the force field of all the other dipoles of the system. The dynamics in this simple model is described essentially by only one parameter eta = 1/3-rho-mu-2/k(B)T, where rho is the number density, mu is the dipole moment, and k(B)T is the Boltzmann constant times the temperature. Extensive Brownian dynamics simulations are carried out to obtain the frequency-dependent dielectric function of this system. The dielectric function becomes progressively non-Debye as the polarity of the system is increased. A comparison between single particle and collective moment-moment time correlation functions is made. It is found that at high polarities, the collective correlation function decays much faster than the single particle function, although a slowly decaying component develops for the collective function. The projection operator technique is used to derive a perturbative equation of motion for the dipole moment at any lattice site-the equation is exact to the second order in eta. This approach is quite successful in explaining the single particle dynamics. The perturbative equation also accounts for the progressively faster initial decay of the collective correlation function as eta increases. However, it fails to describe the slowly decaying component and consequently the observed non-Debye behavior at large eta. We also make a detailed comparison of the simulation results with the self-consistent continuum model of Nee and Zwanzig. It is found that the generalized diffusion equation is reliable if an accurate dielectric friction is used (from simulation, for example). However, the other ingredients of the continuum model are inadequate to describe the dielectric relaxation in this system. In particular, the predicted torque-torque correlation function decays too rapidly compared to the simulation results.