Zirconia has many important phases with Zr coordination varying from six-fold in the orthorhombic phase to eight-fold in the cubic and tetragonal phases. Development of empirical potentials to describe these zirconia phases is an important but long-standing challenge, and it is a bottleneck for theoretical investigation of large zirconia structures. Here, instead of using the standard core–shell model, we developed a new potential for zirconia by combining the long-range Coulomb interaction and bondorder Tersoff model. The bond-order characteristic of the Tersoff potential enables it to be well suited to describe the zirconia phases with different coordination numbers. In particular, the complex monoclinic phase with two inequivalent oxygen atoms, which is difficult to describe with most existing empirical potentials, is well described by this newly developed potential. This potential provides reasonable predictions of most of the static and dynamic properties of various zirconia phases. Besides its clear physical essence, this potential is at least one order of magnitude faster than core–shell based potentials in molecular dynamics simulation. This is because it does not include an ultralight shell that requires an extremely small time step. We also provide potential scripts for the widely used simulation packages GULP and LAMMPS.