On the basis of the semi-empirical rheology model of soft mud under such oscillatory external forces as water waves, a vertical 2-D numerical model is developed for simulating the interaction between surface water waves and underlying soft mud beds. The rheology model used is featured by remarkable hysteresis loops in the relationships of shear strain versus shear stress and the rate of shear strain versus shear stress and strong nonlinearity. A logarithmic grid from the water-mud interface in the vertical direction is employed for numerical treatment, which increases the resolution of the flow in the neighborhood of both sides of the interface, thus contributes to reproducing the characteristics of wave boundary. Model verifications are given through comparisons between the calculated and the measured bed mud mass transport velocities as well as the rate of decay of surface wave heights.