Again: Non-relativistic quantum mechanics IS COVARIANT under Galileo transformations. Otherwise it wouldn't make sense!
Take one free particle as the most simple case. Then the symmetry is realized as
[tex]t'=t, \quad \vec{x}'=\vec{x}-\vec{v} t, \quad \vec{p}'=\vec{p}-m \vec{v}, \psi'(t',\vec{x}')=\exp(\mathrm{i} \Phi) \psi(t,\vec{x})=\exp(\mathrm{i} \Phi) \psi(t',\vec{x}'+\vec{v} t').[/tex]
Now we have to check that we can find a real phase [itex]\Phi(t',\vec{x}')[/itex] such that [itex]\psi'[/itex] fulfills the Schrödinger equation in the new variables, if [itex]\psi(t,\vec{x})[/itex] fulfills it in the old variables. We have
[tex]\mathrm{i} \hbar \partial_{t'} \psi'=\mathrm{i} \hbar \exp(\mathrm{i} \Phi) \left [\mathrm{i} \partial_{t'} \Phi(t',\vec{x}') + \vec{v} \cdot \vec{\nabla} \psi(t,\vec{x})|_{t=t',\vec{x}=\vec{x}'+\vec{v} t'} \right ].[/tex]
Further we find
[tex]\Delta' \psi'(t',\vec{x}') = \left [\mathrm{i} \Delta' \Phi(t',x') \psi(t',\vec{x}'+\vec{v} t') +2 \mathrm{i} \vec{\nabla} ' \Phi(t',\vec{x}') \cdot \vec{\nabla}' \psi(t',\vec{x}'+\vec{v} t') + \Delta ' \psi(t',\vec{x}'+\vec{v} t') + (\vec{\nabla}' \Phi)^2 \psi(t',\vec{x}'+\vec{v} t') \right ]\exp(\mathrm{i} \Phi)[/tex]
Now you demand that
[tex]\mathrm{i} \hbar \partial_{t'} \psi(t',\vec{x}') = -\frac{\hbar^2}{2m} \Delta' \psi(t',\vec{x}'),[/tex]
provided that
[tex]\mathrm{i} \hbar \partial_{t} \psi(t,\vec{x}) = -\frac{\hbar^2}{2m} \Delta \psi(t,\vec{x}).[/tex]
Comparing the coefficients in front of [itex]\vec{\nabla} \psi[/itex] and [itex]\psi[/itex], you get the system of PDEs for the phase factors:
[tex]\vec{\nabla}' \Phi(t',\vec{x}') = -\frac{m}{\hbar} \vec{v}, \quad \partial_{t'} \Phi(t',\vec{x}') = \frac{\hbar}{2m} \mathrm{i} \Delta' \Phi(t',\vec{x}') -\frac{\hbar}{2m} [\vec{\nabla}' \phi(t',\vec{x}')]^2.[/tex]
From the first equation you get
[tex]\Phi(t',\vec{x}')=-\frac{m}{\hbar} \vec{v} \cdot \vec{x} + \Phi_1(t')[/tex]
and from the 2nd
[tex]\dot{\Phi}_1(t')=-\frac{\hbar}{2m} \vec{v}^2 \; \Rightarrow \; \Phi=-\frac{m}{\hbar} \vec{v} \cdot \vec{x}' - \frac{\hbar}{2m} \vec{v}^2 t' + \phi_0.[/tex]
Thus the transformation law for the Schrödinger wave function under Galileo boosts reads
[tex]\psi'(t',x')=\exp\left [\frac{m}{\hbar} \vec{v} \cdot \vec{x}' - \frac{\hbar}{2m} \vec{v}^2 t' + \phi_0 \right ] \psi(t',\vec{x}'+\vec{v} t').[/tex]
This means that up to a phase factor the wave function is a scalar under Galileo boosts and thus quantum theory is covariant under Galileo boosts, because the measurable quantities are all invariant under the Galileo transformation.
The Galileo group is realized as a unitary ray representation. There is no sensible realization of the Galileo group as a unitary transformation. The mass is a non-trivial central charge of the Galileo algebra, and setting it to 0 does not lead to useful quantum dynamics. Massless particles make no sense in non-relativistic physics at all. For a deeper understanding, see
Inönü, E., Wigner, E. P.: Representations of the Galilei group, Il Nuovo Cimento 9(8), 705–718, 1952
http://dx.doi.org/10.1007/BF02782239
Ballentine, Leslie E.: Quantum Mechanics, World Scientific, 1998