In the paper, the inclusive production of heavy quarkonium $eta_Q$ ($Q=b$ or $c$) via $Z$ boson decays within the framework of nonrelativistic QCD (NRQCD) effective theory are studied. The contributions from the leading color-singlet and color-octet Fock states are considered. Total and differential decay widths for the inclusive decays $Z to eta_c+X$ and $Z to eta_b+X$ are presented. It is found that the decay $Zto eta_Q +X$ is dominated by the $^3S_1^{[8]}$ component, so the decays can be inversely adopted to determine the values of the long-distance matrix elements $langle {cal O}^{eta_{c}}(^{3}S_{1}^{[8]})rangle$ and $langle {cal O}^{eta_{b}}(^{3}S_{1}^{[8]})rangle$ respectively. The numerical results show that at an $e^+e^-$ collider running at the $Z$ pole with a high luminosity around $10^{35}{rm cm}^{-2}{rm s}^{-1}$ (a super $Z$ factory), there are about $4.5times 10^7$ $eta_c$ meson events and $6.2times 10^5$ $eta_b$ meson events to be produced per operation year, and the inclusive decays may be used for clarifying some problems on the heavy quarkonium $eta_Q$ and NRQCD.