  We present a very large high-resolution cosmological N-body simulation, the
Millennium-XXL or MXXL, which uses 303 billion particles to represent the
formation of dark matter structures throughout a 4.1Gpc box in a LambdaCDM
cosmology. We create sky maps and identify large samples of galaxy clusters
using surrogates for four different observables: richness estimated from galaxy
surveys, X-ray luminosity, integrated Sunyaev-Zeldovich signal, and lensing
mass. The unprecedented combination of volume and resolution allows us to
explore in detail how these observables scale with each other and with cluster
mass. The scatter correlates between different mass-observable relations
because of common sensitivities to the internal structure, orientation and
environment of clusters, as well as to line-of-sight superposition of
uncorrelated structure. We show that this can account for the apparent
discrepancies uncovered recently between the mean thermal SZ signals measured
for optically and X-ray selected clusters by stacking data from the Planck
satellite. Related systematics can also affect inferences from extreme clusters
detected at high redshift. Our results illustrate that cosmological conclusions
from galaxy cluster surveys depend critically on proper modelling, not only of
the relevant physics, but also of the full distribution of the observables and
of the selection biases induced by cluster identification procedures.
