  We present a new Particle-Mesh cosmological N-body code for accurately
solving the modified Poisson equation of the Quasi Linear formulation of MOND.
We generate initial conditions for the Angus (2009) cosmological model, which
is identical to LCDM except that the cold dark matter is switched for a single
species of thermal sterile neutrinos. We set the initial conditions at z=250
for a (512 Mpc/h)^3 box with 256^3 particles and we evolve them down to z=0. We
clearly demonstrate the necessity of MOND for developing the large scale
structure in a hot dark matter cosmology and contradict the naive expectation
that MOND cannot form galaxy clusters. We find that the correct order of
magnitude of X-ray clusters (with T_X > 4.5 keV) can be formed, but that we
overpredict the number of very rich clusters and seriously underpredict the
number of lower mass clusters. The latter is a shortcoming of the resolution of
our simulations, whereas we suggest that the over production of very rich
clusters might be prevented by incorporating a MOND acceleration constant that
varies with redshift and an expansion history that cannot be described by the
usual Friedmann models. We present evidence that suggests the density profiles
of our simulated clusters are compatible with those of observed X-ray clusters
in MOND. It remains to be seen if the low mass end of the cluster mass function
can be reproduced and if the high densities of dark matter in the central 20
kpc of groups and clusters of galaxies, measured in the MOND framework, can be
achieved. As a last test, we computed the relative velocity between pairs of
halos within 10 Mpc and find that pairs with velocities larger than 3000 km/s
like the bullet cluster, can form without difficulty.
