  We propose a new simple and explicit numerical scheme for time-homogeneous
stochastic differential equations. The scheme is based on sampling increments
at each time step from a skew-symmetric probability distribution, with the
level of skewness determined by the drift and volatility of the underlying
process. We show that as the step-size decreases the scheme converges weakly to
the diffusion of interest. We then consider the problem of simulating from the
limiting distribution of an ergodic diffusion process using the numerical
scheme with a fixed step-size. We establish conditions under which the
numerical scheme converges to equilibrium at a geometric rate, and quantify the
bias between the equilibrium distributions of the scheme and of the true
diffusion process. Notably, our results do not require a global Lipschitz
assumption on the drift, in contrast to those required for the Euler--Maruyama
scheme for long-time simulation at fixed step-sizes. Our weak convergence
result relies on an extension of the theory of Milstein \& Tretyakov to
stochastic differential equations with non-Lipschitz drift, which could also be
of independent interest. We support our theoretical results with numerical
simulations.
