  One of the most challenging open problems in heavy quarkonium physics is the
double charm production in $e^+e^-$ annihilation at B factories. The measured
cross section of $e^+ e^- \to J/\psi + \eta_c$ is much larger than leading
order (LO) theoretical predictions. With the nonrelativistic QCD factorization
formalism, we calculate the next-to-leading order (NLO) QCD correction to this
process. Taking all one-loop self-energy, triangle, box, and pentagon diagrams
into account, and factoring the Coulomb-singular term into the $c\bar c$ bound
state wave function, we get an ultraviolet and infrared finite correction to
the cross section of $e^+e^-\to J/\psi + \eta_c$ at $\sqrt{s} =10.6$ GeV. We
find that the NLO QCD correction can substantially enhance the cross section
with a K factor (the ratio of NLO to LO) of about 1.8-2.1; hence it greatly
reduces the large discrepancy between theory and experiment. With $m_c=1.4{\rm
GeV}$ and $\mu=2m_c$, the NLO cross section is estimated to be 18.9 fb, which
reaches to the lower bound of experiment.
