ABSTRACT As the earliest relics of star formation episodes of the Universe, the most massive galaxies are the key to our understanding of the stellar population, cosmic structure, and supermassive black hole (SMBH) evolution. However, the details of their formation histories remain uncertain. We address these problems by planning a large survey sample of 101 ultramassive galaxies (z ≤ 0.3, |δ + 24°| < 45°, |b| > 8°), including 76 per cent ellipticals, 17 per cent lenticulars, and 7 per cent spirals brighter than MK ≤ −27 mag (stellar mass 2 × 1012 ≲ M⋆ ≲ 5 × 1012 M⊙) with ELT/HARMONI. Our sample comprises diverse galaxy environments ranging from isolated to dense-cluster galaxies. The primary goals of the project are to (1) explore the stellar dynamics inside galaxy nuclei and weigh SMBHs, (2) constrain the black hole scaling relations at the highest mass, and (3) probe the late-time assembly of these most massive galaxies through the stellar population and kinematical gradients. We describe the survey, discuss the distinct demographics and environmental properties of the sample, and simulate their HARMONI Iz-, Iz + J-, and H + K-band observations by combining the inferred stellar-mass models from Pan-STARRS observations, an assumed synthetic spectrum of stars, and SMBHs with masses estimated based on different black hole scaling relations. Our simulations produce excellent state-of-the-art integral field spectrography and stellar kinematics (ΔVrms ≲ 1.5 per cent) in a relatively short exposure time. We use these stellar kinematics in combination with the Jeans anisotropic model to reconstruct the SMBH mass and its error using a Markov chain Monte Carlo simulation. Thus, these simulations and modellings can be benchmarks to evaluate the instrument models and pipelines dedicated to HARMONI to exploit the unprecedented capabilities of ELT.