Based on fast multipole boundary element method (FM-BEM) and mixed variational inequality, a new numerical method named mixed fast multipole boundary element method (MFM-BEM) was presented in this paper for solving three-dimensional elastic-plastic contact problems. Mixed boundary integral equation (MBIE) was the foundation of MFM-BEM and obtained by mixed variational inequality. In order to adapt the requirement of fast multipole method (FMM), Taylor series expansion was used in discrete MBIE. In MFM-BEM the calculation time was significant decreased, the calculation accuracy and continuity was also improved. These merits of MFM-BEM were demonstrated in numerical examples. MFM-BEM has broad application prospects and will take an important role in solving large-scale engineering problems.