A theoretical framework to include the influences of nonbreaking surface waves in ocean general circulation models is established based on Reynolds stresses and fluxes terms derived from surface wave-induced fluctuation. An expression for the wave-induced viscosity and diffusivity as a function of the wave number spectrum is derived for infinite and finite water depths; this derivation allows the coupling of ocean circulation models with a wave number spectrum numerical model. In the case of monochromatic surface wave, the wave-induced viscosity and diffusivity are functions of the Stokes drift. The influence of the wave-induced mixing scheme on global ocean circulation models was tested with the Princeton Ocean Model, indicating significant improvement in upper ocean thermal structure and mixed layer depth compared with mixing obtained by the Mellor-Yamada scheme without the wave influence. For example, the model observation correlation coefficient of the upper 100-m temperature along 35 degrees N increases from 0.68 without wave influence to 0.93 with wave influence. The wave-induced Reynolds stress can reach up to about 5% of the wind stress in high latitudes, and drive 2-3 Sv transport in the global ocean in the form of mesoscale eddies with diameter of 500-1,000 km. The surface wave-induced mixing is more pronounced in middle and high latitudes during the summer in the Northern Hemisphere and in middle latitudes in the Southern Hemisphere.