  The nuclear many-body problem for medium-mass systems is commonly addressed
using wave-function expansion methods that build upon a second-quantized
representation of many-body operators with respect to a chosen computational
basis. While various options for the computational basis are available,
perturbatively constructed natural orbitals recently have been shown to lead to
significant improvement in many-body applications yielding faster model-space
convergence and lower sensitivity to basis set parameters in large-scale
no-core shell model diagonalizations. This work provides a detailed comparison
of single-particle basis sets and a systematic benchmark of natural orbitals in
nonperturbative many-body calculations using the in-medium similarity
renormalization group approach. As a key outcome we find that the construction
of natural orbitals in a large single-particle basis enables for performing the
many-body calculation in a reduced space of much lower dimension, thus offering
significant computational savings in practice that help extend the reach of ab
initio methods towards heavier masses and higher accuracy.
