The contact and non-contact behaviors of material interfaces generate harmonics arising from interactions with large amplitude incident waves. In this phenomenon, which is called contact acoustic nonlinearity (CAN), intermittent impacts by an incident wave cause the interface contact phase to alternate between contact and separation. This paper proposes a method for numerical simulation of CAN using an elastodynamic finite integration technique (EFIT). The accuracy of the EFIT simulation was verified through comparison with an analytical solution. In one-dimensional CAN theory, a wave penetrating at an interface demonstrates a sawtooth waveform indicating that the interface is closing with a constant velocity. Simulation results revealed that the closing velocity is determined by the compressive stress of the material, and experimental measurements on a polymethylmethacrylate specimen revealed the harmonics caused by the sawtooth wave.