  We present a fully differential next-to-next-to-leading order calculation of
the Higgs pair production in association with a $Z$ boson, which is important
for probing the trilinear Higgs self-coupling. The
next-to-next-to-leading-order corrections enhance the next-to-leading order
total cross sections by a factor of $1.2\sim 1.5$, depending on the collider
energy, and change the shape of next-to-leading order kinematic distributions.
We discuss how to determine the trilinear Higgs self-coupling using our
results.
