Shear wave elastography is a promising approach to characterize mechanical properties of soft tissue. By using radiation pressure, it is possible to create deeply localized stress in biological tissue which generates shear waves. We previously studied analytically these induced displacements, supposing that the tissue can be modeled as a linear isotropic solid. Considering a homogeneous solid, we have described the directivities of the compression and shear waves generated by this localized force using the elastodynamic Green's function. If we now consider a heterogeneous solid, the Green's function theory cannot be used. Then we have developed a numerical technique based on a pseudospectral (PS) method, which solves the wave equation with a source term in an isotropic solid. Perfectly matched layers (PML) are used on the boundaries of the numerical grid to avoid reflections. Moreover, the temporal evolution is done with an Adams-Bashforth method. In the case of a point source in a homogeneous elastic medium, numerical results are favorably compared (directivity patterns, displacement shape, ...) with analytical results computed with the Green's function. The shear wave propagation has then been simulated in a 2D plane comprising a small inclusion with a shear elasticity coefficient different from the surrounding tissue. The shear wavefront deformation, as well as the diffraction and reflection of the shear,wave, can be calculated and analyzed.