By performing a 100 ns molecular dynamics simulation of a dipalmitoylphosphatidylcholine lipid bilayer we are able to calculate the full rotational correlation functions of the hydrocarbon chain C-H vectors and determine rotational diffusion of entire lipid molecules with high accuracy. The simulated relaxation is strongly nonexponential already on time scales from 0.1 ps. Fourier transformation of the correlation functions yields data that in the relatively narrow frequency range accessible to H-2 and C-13 NMR experiments are consistent with the reported 1/root omega dependence of spin-lattice relaxation rates. The simulated relaxation dynamics is found to be slightly faster than experimental, which we suggest is explained by the limited accuracy in dihedral potentials of present force fields. By introducing a local frame of reference, the chain motion is separated into local dihedral transitions and overall lipid reorientation. The internal chain isomerization dominates the relaxation and is well-described by power laws. The small molecular reorientation contribution to the decay is exponential with separate time scales for motions of the lipid long axis (D-perpendicular to = 2.9 x 10(7) s(-1)) and spinning rotation around it (D-parallel to = 3.8 x 10(8) s(-1)). The mean square lateral displacement over 100 ns, corrected for the relative motions of the layers, corresponds to a long-time translational diffusion coefficient of D-lat = 1.2 x 10(-7) cm(2) s(-1) at 323 K.