In this paper, incompressible, double-diffusive convection is simulated using finite-difference schemes. The Navier-Stokes equations are expressed in terms of stream function and vorticity. Because of the existence of large velocity, temperature and salinity gradients in boundary layers, a boundary-fitted coordinate system is used to concentrate the grid points near the wall and fit complex boundaries. The finite-difference methods used include the high-order accurate upwind difference scheme. It is shown that the scheme is a good candidate for direct simulations of double-diffusive convection flows. The proposed method is first applied to symmetry breaking and overturning states in thermohaline-driven flows in trapezoid basins. The basic phenomena agree well with those by Dijkstra and Molemaker (1997 J. Fluid Mech. 331 169) and Quon and Ghil (1992 J. Fluid Mech. 245 449), but symmetry breaking and overturning states can occur in an asymmetric geometrical region without perturbations. Then the method is applied to double-diffusive convections in a cavity with opposing horizontal temperature and concentration gradients at large thermal (Rt), solutal (Rs) Rayleigh numbers and Lewis number. There are three straight sides and a sine curve side in the cavity. Basically, numerical results are in agreement with those of Lee and Hyun (1990 Int. J. Heat Mass Transfer 33 1619) qualitatively, but eddies mixing in the top left-hand corner near the curved wall affects the layered structure.