Summary: Coastal flow involves surface wave propagation, current circulation, and seabed evolution, and its prediction remains challenging when they strongly interact with each other, especially during extreme events such as tsunami and storm surge. We propose a fully coupled method to simulate motion of wave-current-seabed systems and associated multiphysics phenomena. The wave action equation, the shallow water equations, and the Exner equation are respectively used for wave, current, and seabed morphology, and the discretization is based on a second-order, flux-limiter, finite difference scheme previously developed for current-seabed systems. The proposed method is tested with analytical solutions, laboratory measurements, and numerical solutions obtained with other schemes. Its advantages are demonstrated in capturing interplay among wave, current, and seabed; it has the capability of first-order upwind schemes to suppress artificial oscillations as well as the accuracy of second-order schemes in resolving flow structures.