  We present a numerical approximation of the solutions of the Euler equations
with a gravitational source term. On the basis of a Suliciu type relaxation
model with two relaxation speeds, we construct an approximate Riemann solver,
which is used in a first order Godunov-type finite volume scheme. This scheme
can preserve both stationary solutions and the low Mach limit to the
corresponding incompressible equations. In addition, we prove that our scheme
preserves the positivity of density and internal energy, that it is entropy
satisfying and also guarantees not to give rise to numerical checkerboard modes
in the incompressible limit. Later we give an extension to second order that
preserves positivity, asymptotic-preserving and well-balancing properties.
Finally, the theoretical properties are investigated in numerical experiments.
