We have developed a two-dimensional P-SV viscoelastic finite-difference modeling technique for complex surface topography and subsurface structures. Realistic modeling of seismic wave propagation in the near surface region is complicated by many factors, such as strong heterogeneity, topographic relief, and large attenuation. In order to account for these complications, we use an O(2,4) accurate viscoelastic velocity-stress staggered-grid finite-difference scheme. The implementation includes an irregular free surface condition for topographic relief and a discontinuous grid technique in the shallow parts of the model. Several methods of free surface condition are bench marked, and an accurate and simple condition is proposed. In the proposed free surface condition, stresses are calculated so that the normal stresses perpendicular to the boundary and shear stresses on the free surface are zero. The calculation of particle velocities at the free surface does not involve any specific calculations, and the particle velocities are set to zero above the free surface. A discontinuous-grid method is introduced, where we use a 3 times finer grid in the near surface or low velocity region compared to the rest of the model. In order to reduce instability, we apply averaging or weighting to the replacement of the coarse-grid components within the fine grid field. The method allows us to avoid any limitation of the shape of the grid-spacing boundary. Numerical tests indicate that approximately 10 grid points per shortest wavelength, counted in coarse-grid spacing, with the discontinuous grid method results in accurate calculations as long as a small number of time steps is concerned.