  Density functional theory (DFT) has become a standard tool for the study of
point defects in materials. However, finding the most stable defective
structures remains a very challenging task as it involves the solution of a
multimodal optimization problem with a high-dimensional objective function.
Hitherto, the approaches most commonly used to tackle this problem have been
mostly empirical, heuristic and/or based on domain knowledge. In this
contribution, we describe an approach for exploring the potential energy
surface based on the covariance matrix adaption evolution strategy (CMA-ES) and
supervised and unsupervised machine learning models. We show how the original
CMA-ES can be modified to suit the specific problem of DFT studies of point
defects in the dilute limit. The resulting algorithm depends only on a limited
set of physically interpretable hyperparameters. The approach offers a robust
and systematic way for finding low-energy configurations of point defects in
solids. We demonstrate the applicability and moderate computational cost on the
intrinsic defects in silicon. We also apply the methodology to the neutral
oxygen vacancy oxygen vacancy in TiO$_2$ anatase and reproduce the known defect
structures. Furthermore, a new defect structure, stable at the level of hybrid
density functional theory and characterized by a delocalized electronic
structure, is found for this system.
