Abstract: A model of gear contact has been established based on Hertz theory. Local coordinates are constructed on different meshing points because the curvature of gear profile is not constant, and then distributions of residual stresses and strains are constructed. A numerical approach for elastic-plastic analysis and shakedown is presented and its application to gear contact is demonstrated. The numerical approach consists of an operator split technique, which transforms the elastic-plastic problem into a purely elastic problem and a residual problem with prescribed eigenstrains. The eigenstrains are determined using an incremental projection method. The proposed method is applied to compute contact stresses and contact residual stresses of meshing gear teeth with standard and modified profile. The comparison with experimental results demonstrates the validity of the proposed method.