  In this paper, we propose and analyse a numerical method to solve 2D
Dirichlet time-harmonic elastic wave equations. The procedure is based on the
decoupling of the elastic vector field into scalar Pressure ($P$-) and Shear
($S$-) waves via a suitable Helmholtz-Hodge decomposition. For the
approximation of the two scalar potentials we apply a virtual element method
associated with different mesh sizes and degrees of accuracy. We provide for
the stability of the method and a convergence error estimate in the $L^2$-norm
for the displacement field, in which the contributions to the error associated
with the $P$- and $S$- waves are separated. In contrast to standard approaches
that are directly applied to the vector formulation, this procedure allows for
keeping track of the two different wave numbers, that depend on the $P$- and
$S$- speeds of propagation and, therefore, for using a high-order method for
the approximation of the wave associated with the higher wave number. Some
numerical tests, validating the theoretical results and showing the good
performance of the proposed approach, are presented.
