  We enhance the treatment of crystallization for models of white dwarfs (WDs)
in the stellar evolution software MESA by implementing carbon-oxygen (C/O)
phase separation. The phase separation process during crystallization leads to
transport of oxygen toward the center of WDs, resulting in a more compact
structure that liberates gravitational energy as additional heating that
modestly slows WD cooling timescales. We quantify this cooling delay in MESA
C/O WD models over the mass range 0.5-1.0 $M_\odot$, finding delays of 0.5-0.8
Gyr for typical C/O interior profiles. MESA WD cooling timescales including
this effect are generally comparable to other WD evolution models that make
similar assumptions about input physics. When considering phase separation
alongside $^{22}$Ne sedimentation, however, we find that both MESA and BaSTI WD
cooling models predict a more modest sedimentation delay than the latest LPCODE
models, and this may therefore require a re-evaluation of previously proposed
solutions to some WD cooling anomalies that were based on LPCODE models of
$^{22}$Ne sedimentation. Our implementation of C/O phase separation in the
open-source stellar evolution software MESA provides an important tool for
building realistic grids of WD cooling models, as well as a framework for
expanding on our implementation to explore additional physical processes
related to phase transitions and associated fluid motions in WD interiors.
